#include "mpi.h"#include "stdlib.h"#include "emdata.h"#include "alignoptions.h"Include dependency graph for utilcomm_Cart.h:

This graph shows which files directly or indirectly include this file:

Go to the source code of this file.
Functions | |
| int | ReadStackandDist_Cart (MPI_Comm comm, EMData ***images2D, char *stackfname, int *nloc) |
| int | CleanStack_Cart (MPI_Comm comm, EMData **image_stack, int nloc, int ri, Vec3i volsize, Vec3i origin) |
| int | ReadAngTrandDist_Cart (MPI_Comm comm_2d, MPI_Comm comm_row, int *dim, float *angleshift, char *angfname, int nloc) |
| int | setpart_gc1 (MPI_Comm comm_2d, int nangs, int *psize, int *nbase) |
| int | setpart_gr1 (MPI_Comm comm_2d, int nnz, int *nnzpart, int *nnzbase) |
| int | sphpart (MPI_Comm comm_2d, int nrays, int *ptrs, int *nnzbase, int *ptrstart) |
| int | getcb2sph (Vec3i volsize, int ri, Vec3i origin, int nnz0, int *ptrs, int *cord) |
| int | fwdpj3_Cart (Vec3i volsize, int nraysloc, int nnzloc, float *dm, Vec3i origin, int ri, int *ptrs, int *cord, int myptrstart, float *x, float *y) |
| int | bckpj3_Cart (Vec3i volsize, int nraysloc, int nnzloc, float *dm, Vec3i origin, int ri, int *ptrs, int *cord, int myptrstart, float *x, float *y) |
|
||||||||||||||||||||||||||||||||||||||||||||||||
|
Definition at line 237 of file project3d_Cart.cpp. References cord, dm, ifix(), myptrstart, nx, ny, ptrs, status, EMAN::Vec3i, x, and y. 00240 : 00241 volsize - size of 3D cube volume 00242 nraysloc - local number of rays 00243 nnzloc - local number of voxels 00244 dm - projection matrix 00245 origin, ri - origin and radius of sphere 00246 ptrs, cord - pointers and coordinates for first ray 00247 myptrstart - determines starting rays 00248 x - 2D image 00249 Output: 00250 y - portion of the 3D volume 00251 */ 00252 { 00253 int i, j, iqx,iqy, xc, yc, zc, jj; 00254 float xb, yb, dx, dy, dx1m, dy1m, dxdy; 00255 int status = 0; 00256 00257 int xcent = origin[0]; 00258 int ycent = origin[1]; 00259 int zcent = origin[2]; 00260 00261 int nx = volsize[0]; 00262 int ny = volsize[1]; 00263 00264 // Phi: adding the shift parameters that get passed in as the last two entries of dm 00265 float sx, sy; 00266 00267 sx = dm(7); 00268 sy = dm(8); 00269 00270 if ( nx > 2*ri) { 00271 for (i = 1; i <= nraysloc; i++) { 00272 zc = cord(1,myptrstart+i) - zcent; 00273 yc = cord(2,myptrstart+i) - ycent; 00274 xc = cord(3,myptrstart+i) - xcent; 00275 00276 xb = zc*dm(1)+yc*dm(2)+xc*dm(3) + xcent + sx; 00277 yb = zc*dm(4)+yc*dm(5)+xc*dm(6) + ycent + sy; 00278 00279 for (j = ptrs(myptrstart+i); j <ptrs(myptrstart+i+1); j++) { 00280 jj = j-ptrs[myptrstart]+1; 00281 00282 // jj is the index on the local volume 00283 iqx = ifix(xb); 00284 iqy = ifix(yb); 00285 00286 dx = xb - iqx; 00287 dy = yb - iqy; 00288 dx1m = 1.0 - dx; 00289 dy1m = 1.0 - dy; 00290 dxdy = dx*dy; 00291 /* 00292 c y(j) = y(j) + dx1m*dy1m*x(iqx , iqy) 00293 c & + dx1m*dy *x(iqx , iqy+1) 00294 c & + dx *dy1m*x(iqx+1, iqy) 00295 c & + dx *dy *x(iqx+1, iqy+1) 00296 c 00297 c --- faster version of the above commented out 00298 c code (derived by summing the following table 00299 c of coefficients along the colunms) --- 00300 c 00301 c 1 dx dy dxdy 00302 c ------ -------- -------- ------- 00303 c x(i,j) -x(i,j) -x(i,j) x(i,j) 00304 c x(i,j+1) -x(i,j+1) 00305 c x(i+1,j) -x(i+1,j) 00306 c x(i+1,j+1) 00307 c 00308 */ 00309 // Phi: add index checking, now that shifts are being used 00310 if ( iqx <= nx && iqy <= ny && iqx >= 1 && iqy >= 1 ) { 00311 y(jj) += x(iqx,iqy); 00312 if ( iqx + 1 <= nx && iqx + 1 >= 1 ) { 00313 y(jj) += dx*(-x(iqx,iqy)+x(iqx+1,iqy)); 00314 } 00315 if ( iqy + 1 <= ny && iqy + 1 >= 1 ) { 00316 y(jj) += dy*(-x(iqx,iqy)+x(iqx,iqy+1)); 00317 } 00318 if ( iqx + 1 <= nx && iqy + 1 <= ny && iqx + 1 >= 1 && iqy + 1 >= 1 ) { 00319 y(jj) += dxdy*( x(iqx,iqy) - x(iqx,iqy+1) -x(iqx+1,iqy) + x(iqx+1,iqy+1) ); 00320 } 00321 } 00322 00323 // y(j) += x(iqx,iqy) 00324 // + dx*(-x(iqx,iqy)+x(iqx+1,iqy)) 00325 // + dy*(-x(iqx,iqy)+x(iqx,iqy+1)) 00326 // + dxdy*( x(iqx,iqy) - x(iqx,iqy+1) 00327 // -x(iqx+1,iqy) + x(iqx+1,iqy+1) ); 00328 00329 xb += dm(1); 00330 yb += dm(4); 00331 } // end for j 00332 } // end for i 00333 } 00334 else { 00335 fprintf(stderr, "bckpj3: nx must be greater than 2*ri\n"); 00336 } 00337 00338 return status; 00339 }
|
|
||||||||||||||||||||||||||||
|
Definition at line 275 of file utilcomm_Cart.cpp. References asta2(), EMDeleteArray(), EMAN::EMData::get_data(), imgdata, nx, ny, rhs, and EMAN::Vec3i. Referenced by main(). 00276 {
00277 int nx = volsize[0];
00278 int ny = volsize[1];
00279
00280 float * rhs = new float[nx*ny];
00281 float * imgdata;
00282
00283 // Calculate average "background" from all pixels strictly outside of radius
00284 double aba, abaloc; // average background
00285 aba = 0.0;
00286 abaloc = 0.0;
00287 int klp, klploc; // number of pixels in background
00288 klp = 0;
00289 klploc = 0;
00290
00291 // calculate avg background in parallel
00292 for ( int i = 0 ; i < nloc ; ++i ) {
00293 imgdata = image_stack[i]->get_data();
00294 asta2(imgdata, nx, ny, ri, &abaloc, &klploc);
00295 }
00296
00297 // All reduce using the column-based sub-topology
00298
00299 MPI_Allreduce(&abaloc, &aba, 1, MPI_DOUBLE, MPI_SUM, comm_col);
00300 MPI_Allreduce(&klploc, &klp, 1, MPI_INT, MPI_SUM, comm_col);
00301
00302 aba /= klp;
00303
00304 // subtract off the average background from pixels weakly inside of radius
00305 int x_summand, y_summand;
00306 int r_squared = ri * ri;
00307 for ( int i = 0 ; i < nloc ; ++i ) {
00308 imgdata = image_stack[i]->get_data();
00309 for ( int j = 0 ; j < nx ; ++j) {
00310 x_summand = (j - origin[0]) * (j - origin[0]);
00311 for ( int k = 0 ; k < ny ; ++k ) {
00312 y_summand = (k - origin[1]) * (k - origin[1]);
00313 if ( x_summand + y_summand <= r_squared) {
00314 imgdata[j*ny + k] -= aba;
00315 }
00316 }
00317 }
00318 }
00319 EMDeleteArray(rhs);
00320 return 0;
00321 }
|
|
||||||||||||||||||||||||||||||||||||||||||||||||
|
Definition at line 147 of file project3d_Cart.cpp. References cord, dm, ifix(), myptrstart, nx, ny, ptrs, status, EMAN::Vec3i, x, and y. 00150 :
00151 volsize - size (nx,ny,nz) of 3D cube volume
00152 nraysloc - local number of rays within the compact spherical representation
00153 nnzloc - local number of voxels
00154 dm - an array of size 9 storing transformation
00155 origin, ri - origin and radius of sphere
00156 ptrs, cord - pointers and coordinates for first ray
00157 myptrstart - determines starting rays
00158 x - portion of the 3D input volume
00159 Output:
00160 y - 2D output image
00161 */
00162 {
00163 int iqx, iqy, i, j, xc, yc, zc, jj;
00164 float ct, dipx, dipy, dipx1m, dipy1m, xb, yb, dm1, dm4;
00165 int status = 0;
00166
00167 // Phi: adding the shift parameters that get passed in as the last two entries of dm
00168 float sx, sy;
00169
00170 sx = dm(7);
00171 sy = dm(8);
00172
00173 int xcent = origin[0];
00174 int ycent = origin[1];
00175 int zcent = origin[2];
00176
00177 int nx = volsize[0];
00178 int ny = volsize[1];
00179
00180 dm1 = dm(1);
00181 dm4 = dm(4);
00182
00183 if ( nx > 2*ri ) {
00184 for (i = 1; i <= nraysloc; i++) {
00185
00186 zc = cord(1,myptrstart+i)-zcent;
00187 yc = cord(2,myptrstart+i)-ycent;
00188 xc = cord(3,myptrstart+i)-xcent;
00189 xb = zc* dm(1) +yc* dm(2) +xc* dm(3) + xcent + sx;
00190 yb = zc* dm(4) +yc* dm(5) +xc* dm(6) + ycent + sy;
00191
00192 for (j = ptrs(myptrstart+i); j< ptrs(myptrstart+i+1); j++) {
00193 jj = j-ptrs[myptrstart]+1;
00194 iqx = ifix(xb);
00195 iqy = ifix(yb);
00196
00197 // jj is the index on the local volume
00198 ct = x(jj);
00199
00200 // dipx = xb - (float)(iqx);
00201 // dipy = (yb - (float)(iqy)) * ct;
00202 dipx = xb - iqx;
00203 dipy = (yb - iqy) * ct;
00204
00205 dipy1m = ct - dipy;
00206 dipx1m = 1.0 - dipx;
00207
00208 if (iqx <= nx && iqy <= ny && iqx >= 1 && iqy >= 1)
00209 // y(iqx ,iqy) = y(iqx ,iqy) + dipx1m*dipy1m;
00210 y(iqx ,iqy) += dipx1m*dipy1m;
00211 if (iqx + 1 <= nx && iqy <= ny && iqx >= 0 && iqy >= 1)
00212 // y(iqx+1,iqy) = y(iqx+1,iqy) + dipx*dipy1m;
00213 y(iqx+1,iqy) += dipx*dipy1m;
00214 if (iqx + 1 <= nx && iqy + 1 <= ny && iqx >= 0 && iqy >= 0)
00215 // y(iqx+1,iqy+1) = y(iqx+1,iqy+1) + dipx*dipy;
00216 y(iqx+1,iqy+1) += dipx*dipy;
00217 if (iqx <= nx && iqy + 1 <= ny && iqx >= 1 && iqy >= 0)
00218 // y(iqx ,iqy+1) = y(iqx ,iqy+1) + dipx1m*dipy;
00219 y(iqx ,iqy+1) += dipx1m*dipy;
00220 xb += dm1;
00221 yb += dm4;
00222 }
00223 }
00224 }
00225 else {
00226 fprintf(stderr, " nx must be greater than 2*ri\n");
00227 exit(1);
00228 }
00229 return status;
00230 }
#undef x
|
|
||||||||||||||||||||||||||||
|
Definition at line 15 of file project3d_Cart.cpp. References cord, nnz, nrays, nx, ny, ptrs, status, and EMAN::Vec3i. 00018 : 00019 volsize - dimensions of 3D cube 00020 ri - radius of sphere 00021 origin - origin of sphere 00022 nnz0 - number of non-zeros 00023 00024 Output: 00025 ptrs - vector containing pointers to the beginning of each ray 00026 cord - vector containing the x,y,z coordinates of the beginning of each ray 00027 */ 00028 00029 { 00030 int xs, ys, zs, xx, yy, zz, rs, r2; 00031 int ix, iy, iz, jnz, nnz, nrays; 00032 int ftm = 0, status = 0; 00033 00034 int xcent = (int)origin[0]; 00035 int ycent = (int)origin[1]; 00036 int zcent = (int)origin[2]; 00037 00038 int nx = (int)volsize[0]; 00039 int ny = (int)volsize[1]; 00040 int nz = (int)volsize[2]; 00041 00042 r2 = ri*ri; 00043 nnz = 0; 00044 nrays = 0; 00045 ptrs(1) = 1; 00046 00047 for (ix = 1; ix <= nx; ix++) { 00048 xs = ix-xcent; 00049 xx = xs*xs; 00050 for ( iy = 1; iy <= ny; iy++ ) { 00051 ys = iy-ycent; 00052 yy = ys*ys; 00053 jnz = 0; 00054 00055 ftm = 1; 00056 // not the most efficient implementation 00057 for (iz = 1; iz <= nz; iz++) { 00058 zs = iz-zcent; 00059 zz = zs*zs; 00060 rs = xx + yy + zz; 00061 if (rs <= r2) { 00062 jnz++; 00063 nnz++; 00064 // sphere(nnz) = cube(iz, iy, ix); 00065 00066 // record the coordinates of the first nonzero === 00067 if (ftm) { 00068 nrays++; 00069 cord(1,nrays) = iz; 00070 cord(2,nrays) = iy; 00071 cord(3,nrays) = ix; 00072 ftm = 0; 00073 } 00074 } 00075 } // end for (iz..) 00076 if (jnz > 0) { 00077 ptrs(nrays+1) = ptrs(nrays) + jnz; 00078 } // endif (jnz) 00079 } // end for iy 00080 } // end for ix 00081 if (nnz != nnz0) status = -1; 00082 return status; 00083 }
|
|
||||||||||||||||||||||||||||
|
Definition at line 187 of file utilcomm_Cart.cpp. References COL, EMDeleteArray(), ierr, and ROW. Referenced by main(). 00188 {
00189 int ierr = 0;
00190 int ROW = 0, COL = 1;
00191 int my2dpid, mycoords[2];
00192 int srpid, srcoords[2], keep_dims[2];
00193 MPI_Status mpistatus;
00194 FILE *fp = NULL;
00195
00196 int nimgs=0;
00197
00198 MPI_Comm_rank(comm_2d, &my2dpid); //Get my pid in the new 2D topology
00199 MPI_Cart_coords(comm_2d, my2dpid, 2, mycoords); // Get my coordinates
00200
00201 float * iobuffer = new float[5*nloc];
00202 if (!iobuffer) {
00203 fprintf(stderr,"failed to allocate buffer to read angles shifts\n");
00204 ierr = -1;
00205 goto EXIT;
00206 }
00207
00208 if (mycoords[COL] == 0 && mycoords[ROW] == 0) { //I am Proc (0,0)
00209 fp = fopen(angfname,"r");
00210 if (!fp) ierr = 1;
00211 }
00212 MPI_Bcast(&ierr, 1, MPI_INT, 0, comm_2d);
00213
00214 if ( ierr ) {
00215 if (mycoords[COL] == 0 && mycoords[ROW] == 0)
00216 fprintf(stderr,"failed to open %s\n", angfname);
00217 ierr = MPI_Finalize();
00218 goto EXIT;
00219 }
00220 else {
00221 if (mycoords[COL] == 0 && mycoords[ROW] == 0) { //I am Proc (0,0)
00222 for (int iproc = 0; iproc < dims[ROW]; iproc++) {
00223 // figure out the number of images assigned to processor (iproc,0)
00224 if (iproc > 0) {
00225 srcoords[COL] = 0;
00226 srcoords[ROW] = iproc;
00227 MPI_Cart_rank(comm_2d, srcoords, &srpid);
00228
00229 MPI_Recv(&nimgs, 1, MPI_INT, srpid, srpid, comm_2d, &mpistatus);
00230
00231 // Read the next nimgs set of angles and shifts
00232 for (int i = 0; i < nimgs; i++) {
00233 fscanf(fp,"%f %f %f %f %f",
00234 &iobuffer[5*i+0],
00235 &iobuffer[5*i+1],
00236 &iobuffer[5*i+2],
00237 &iobuffer[5*i+3],
00238 &iobuffer[5*i+4]);
00239 }
00240 MPI_Send(iobuffer,5*nimgs,MPI_FLOAT,srpid,srpid,comm_2d);
00241 }
00242 else {
00243 for (int i = 0; i < nloc; i++) {
00244 fscanf(fp,"%f %f %f %f %f",
00245 &angleshift[5*i+0],
00246 &angleshift[5*i+1],
00247 &angleshift[5*i+2],
00248 &angleshift[5*i+3],
00249 &angleshift[5*i+4]);
00250 }
00251 }
00252 }
00253 fclose(fp);
00254 }
00255 else if (mycoords[COL] == 0 && mycoords[ROW] != 0) { //I am in the first column
00256 // send image count to the master processor (mypid = 0)
00257
00258 MPI_Send(&nloc, 1, MPI_INT, 0, my2dpid, comm_2d);
00259 // Receive angleshifts
00260 MPI_Recv(angleshift, 5*nloc, MPI_FLOAT, 0, my2dpid, comm_2d, &mpistatus);
00261 }
00262 }
00263
00264 // Now have all the processors in group g_c_0 broadcast the angles along the row communicator
00265 srcoords[ROW] = 0;
00266 MPI_Cart_rank(comm_row, srcoords, &srpid);
00267 MPI_Bcast(angleshift, 5*nloc, MPI_FLOAT, srpid, comm_row);
00268
00269 EMDeleteArray(iobuffer);
00270
00271 EXIT:
00272 return ierr;
00273 }
|
|
||||||||||||||||||||
|
Definition at line 11 of file utilcomm_Cart.cpp. References COL, EMDeleteArray(), EMDeletePtr(), EMAN::EMData::get_attr(), EMAN::EMData::get_data(), EMAN::EMUtil::get_image_count(), EMAN::EMData::get_xsize(), EMAN::EMData::get_ysize(), ierr, img_ptr, imgdata, nx, ny, EMAN::EMData::read_image(), ROW, and setpart_gc1(). Referenced by main(). 00012 {
00013 int ncpus, ierr, mpierr=0;
00014 MPI_Status mpistatus;
00015 EMUtil *my_util;
00016 int nima;
00017 FILE *fp=NULL;
00018 int ROW = 0, COL = 1;
00019 int my2dpid, mycoords[2], dims[2], periods[2];
00020 int srpid, srcoords[2]; // Send/receive info
00021
00022 // Get dims and my coordinates
00023 MPI_Cart_get(comm_2d, 2, dims, periods, mycoords);
00024 MPI_Comm_rank(comm_2d, &my2dpid); //Get my pid in the new 2D topology
00025 ncpus = dims[ROW]*dims[COL];
00026
00027 ierr = 0;
00028 if (mycoords[ROW] == 0 && mycoords[COL] == 0) {
00029 fp = fopen(stackfname,"r");
00030 if (!fp) {
00031 ierr = 1;
00032 printf("failed to open %s\n", stackfname);
00033 }
00034 else {
00035 fclose(fp);
00036 nima = my_util->get_image_count(stackfname);
00037 /**********************************************/
00038 // nima = 84;
00039 /****************************/
00040 }
00041 }
00042 mpierr = MPI_Bcast(&ierr, 1, MPI_INT, 0, comm_2d);
00043 if (ierr == 0) {
00044 int *psize = new int[dims[ROW]];
00045 int *nbase = new int[dims[ROW]];
00046
00047 mpierr = MPI_Bcast(&nima, 1, MPI_INT, 0, comm_2d);
00048
00049 *nloc = setpart_gc1(comm_2d, nima, psize, nbase);
00050 *images2D = new EMData*[*nloc]; // NB!: whoever calls ReadStackandDist must delete this!
00051
00052 EMData *img_ptr;
00053 int img_index;
00054 float *imgdata;
00055
00056 // read the first image to get size
00057 img_ptr = new EMData();
00058 img_ptr->read_image(stackfname, 0);
00059 int nx = img_ptr->get_xsize();
00060 int ny = img_ptr->get_ysize();
00061
00062 float s2x, s2y;
00063
00064 if (mycoords[COL] == 0 && mycoords[ROW] == 0) {
00065 printf("Master node reading and distributing %d images along the first column of processors\n", nima);
00066 for ( int ip = 0 ; ip < dims[ROW] ; ++ip ) {
00067 for ( int i = 0 ; i < psize[ip] ; ++i ) {
00068 img_index = nbase[ip] + i;
00069
00070 if (ip != 0) {
00071 img_ptr->read_image(stackfname, img_index);
00072 // get a pointer to the image's data
00073 imgdata = img_ptr->get_data();
00074 // find the x/y shift values if it has them, otherwise set them to 0.0
00075 try {
00076 s2x = (*images2D)[i]->get_attr("s2x");
00077 } catch ( std::exception& e ) {
00078 s2x = 0.0;
00079 }
00080 try {
00081 s2y = (*images2D)[i]->get_attr("s2y");
00082 } catch ( std::exception& e ) {
00083 s2y = 0.0;
00084 }
00085 // send these to processor (ip,0)
00086 srcoords[COL] = 0;
00087 srcoords[ROW] = ip;
00088 MPI_Cart_rank(comm_2d, srcoords, &srpid);
00089
00090 MPI_Send(imgdata, nx*ny, MPI_FLOAT, srpid, srpid, comm_2d);
00091 MPI_Send(&s2x, 1, MPI_FLOAT, srpid, srpid, comm_2d);
00092 MPI_Send(&s2y, 1, MPI_FLOAT, srpid, srpid, comm_2d);
00093 } else { // ip == 0
00094 (*images2D)[i] = new EMData();
00095 (*images2D)[i]->read_image(stackfname, img_index);
00096 try {
00097 s2x = (*images2D)[i]->get_attr("s2x");
00098 } catch ( std::exception& e ) {
00099 (*images2D)[i]->set_attr("s2x",0.0);
00100 }
00101 try {
00102 s2y = (*images2D)[i]->get_attr("s2y");
00103 } catch ( std::exception& e ) {
00104 (*images2D)[i]->set_attr("s2y",0.0);
00105 }
00106 }
00107 }
00108 //printf("finished reading data for processor %d\n", ip);
00109 }
00110
00111 } else if (mycoords[COL] == 0 && mycoords[ROW] != 0) {
00112 for ( int i = 0 ; i < psize[mycoords[ROW]] ; ++i ) {
00113 (*images2D)[i] = new EMData();
00114 (*images2D)[i]->set_size(nx, ny, 1);
00115 imgdata = (*images2D)[i]->get_data();
00116
00117 MPI_Recv(imgdata, nx*ny, MPI_FLOAT, 0, my2dpid, comm_2d, &mpistatus);
00118 MPI_Recv(&s2x, 1, MPI_FLOAT, 0, my2dpid, comm_2d, &mpistatus);
00119 MPI_Recv(&s2y, 1, MPI_FLOAT, 0, my2dpid, comm_2d, &mpistatus);
00120 (*images2D)[i]->set_attr("s2x",s2x);
00121 (*images2D)[i]->set_attr("s2y",s2y);
00122 }
00123 //printf("received %d images for processor (%d,%d)\n", *nloc, mycoords[ROW], mycoords[COL]);
00124 }
00125
00126 // Now have all the processors in group g_c_0 broadcast the images and angles along the row communicator
00127 if (mycoords[COL] == 0 && mycoords[ROW] == 0)
00128 printf("First column of processors distributing images to other columns..\n");
00129
00130 if (mycoords[COL] == 0) {
00131 for (int img_index = 0; img_index < *nloc; img_index++){
00132 imgdata = (*images2D)[img_index]->get_data();
00133
00134 // find the x/y shift values if it has them, otherwise set them to 0.0
00135 try {
00136 s2x = (*images2D)[img_index]->get_attr("s2x");
00137 } catch ( std::exception& e ) {
00138 s2x = 0.0;
00139 }
00140 try {
00141 s2y = (*images2D)[img_index]->get_attr("s2y");
00142 } catch ( std::exception& e ) {
00143 s2y = 0.0;
00144 }
00145
00146 for(int ip=1; ip< dims[COL]; ip++){
00147 // printf("Proc (%d, %d) sending image %d of %d to Proc (%d,%d) \n", mycoords[ROW], mycoords[COL], img_index, *nloc, mycoords[ROW], ip);
00148 srcoords[ROW] = mycoords[ROW];
00149 srcoords[COL] = ip;
00150 MPI_Cart_rank(comm_2d, srcoords, &srpid);
00151
00152 MPI_Send(imgdata, nx*ny, MPI_FLOAT, srpid, img_index, comm_2d);
00153 MPI_Send(&s2x, 1, MPI_FLOAT, srpid, img_index, comm_2d);
00154 MPI_Send(&s2y, 1, MPI_FLOAT, srpid, img_index, comm_2d);
00155 }
00156 }
00157 }
00158 if (mycoords[COL] != 0) {
00159 // printf("Proc (%d, %d) receiving...", mycoords[ROW], mycoords[COL]);
00160 for(int img_index = 0; img_index < *nloc; img_index++){
00161 (*images2D)[img_index] = new EMData();
00162 (*images2D)[img_index]->set_size(nx, ny, 1);
00163 imgdata = (*images2D)[img_index]->get_data();
00164
00165 srcoords[ROW] = mycoords[ROW];
00166 srcoords[COL] = 0;
00167 MPI_Cart_rank(comm_2d, srcoords, &srpid);
00168
00169 MPI_Recv(imgdata, nx*ny, MPI_FLOAT, srpid, img_index, comm_2d, &mpistatus);
00170 MPI_Recv(&s2x, 1, MPI_FLOAT, srpid, img_index, comm_2d, &mpistatus);
00171 MPI_Recv(&s2y, 1, MPI_FLOAT, srpid, img_index, comm_2d, &mpistatus);
00172
00173 (*images2D)[img_index]->set_attr("s2x",s2x);
00174 (*images2D)[img_index]->set_attr("s2y",s2y);
00175 }
00176 }
00177
00178 ierr = mpierr;
00179
00180 EMDeletePtr(img_ptr);
00181 EMDeleteArray(psize);
00182 EMDeleteArray(nbase);
00183 }
00184 return ierr;
00185 }
|
|
||||||||||||||||||||
|
Definition at line 323 of file utilcomm_Cart.cpp. 00326 : 00327 comm_2d - Cartesian communicator 00328 nangs - number of 2D images 00329 Output: 00330 psize - vector containing the partition distribution 00331 nbase - vector containing the base of the partitions 00332 nangsloc - number of local angles 00333 */ 00334 { 00335 int ROW = 0, COL = 1; 00336 int dims[2], periods[2], mycoords[2]; 00337 int nangsloc, nrem; 00338 00339 // Get information associated with comm_2d 00340 MPI_Cart_get(comm_2d, 2, dims, periods, mycoords); 00341 00342 nangsloc = nangs/dims[ROW]; 00343 nrem = nangs - dims[ROW]*nangsloc; 00344 if (mycoords[ROW] < nrem) nangsloc++; 00345 00346 for (int i = 0; i < dims[ROW]; i++) { 00347 psize[i] = nangs/dims[ROW]; 00348 if (i < nrem) psize[i]++; 00349 } 00350 00351 for (int i = 0; i < dims[ROW]; i++) { 00352 if(i==0) 00353 nbase[i] = 0; 00354 else 00355 nbase[i] = nbase[i-1] + psize[i-1]; 00356 } 00357 00358 //printf("I am [%d,%d], nloc = %d, psize = [%d, %d], nbase = [%d, %d]", mycoords[ROW], mycoords[COL],nangsloc,psize[0],psize[1], nbase[0], nbase[1]); 00359 00360 return nangsloc; 00361 }
|
|
||||||||||||||||||||
|
Definition at line 364 of file utilcomm_Cart.cpp. 00367 : 00368 comm_2d - Cartesian communicator 00369 nnz - total number of non-zeros in the volume 00370 Output: 00371 nnzpart - vector containing the partition distribution 00372 nnzbase - vector containing the base of the partitions 00373 nnzloc - initial local nnz 00374 */ 00375 { 00376 int ROW = 0, COL = 1; 00377 int dims[2], periods[2], mycoords[2]; 00378 int nnzloc, nrem; 00379 00380 // Get information associated with comm_2d 00381 MPI_Cart_get(comm_2d, 2, dims, periods, mycoords); 00382 00383 nnzloc = nnz/dims[COL]; 00384 nrem = nnz - dims[COL]*nnzloc; 00385 if (mycoords[COL] < nrem) nnzloc++; 00386 00387 for (int i = 0; i < dims[COL]; i++) { 00388 nnzpart[i] = nnz/dims[COL]; 00389 if (i < nrem) nnzpart[i]++; 00390 } 00391 00392 nnzbase[0] = 0; 00393 for (int i = 1; i < dims[COL]; i++) { 00394 nnzbase[i] = nnzbase[i-1] + nnzpart[i-1]; 00395 } 00396 nnzbase[dims[COL]] = nnz; 00397 return nnzloc; 00398 }
|
|
||||||||||||||||||||||||
|
Definition at line 85 of file project3d_Cart.cpp. 00088 : 00089 comm_2d - Cartesian communicator 00090 nrays - total number of rays 00091 ptrs - vector containing pointers 00092 nnzbase - ideal volume partition of nnz 00093 Output: 00094 ptrstart - vector containing all the starting pointers for each column processor group 00095 nraysloc - actual number of local rays 00096 */ 00097 { 00098 int ROW = 0, COL = 1; 00099 int dims[2], periods[2], mycoords[2]; 00100 int nraysloc = 0; 00101 00102 // Get information associated with comm_2d 00103 MPI_Cart_get(comm_2d, 2, dims, periods, mycoords); 00104 00105 int count = 1; 00106 if (mycoords[COL] == 0){ // First column starts out with the first ray 00107 nraysloc++; 00108 } 00109 ptrstart[0] = 0; 00110 00111 for(int i=1; i<nrays ; i++){ 00112 if (ptrs[i]-1 <= nnzbase[count] && ptrs[i+1]-1 >= nnzbase[count]){ 00113 //nnzbase is between or equal to ptrs 00114 00115 if (nnzbase[count]-(ptrs[i]-1)>= ptrs[i+1]-1-nnzbase[count]){ 00116 if(mycoords[COL] == count-1){ // ray belongs to count-1 00117 nraysloc++; 00118 } 00119 ptrstart[count] = i+1; 00120 count++; 00121 00122 } else { //nnzbase[count]-(ptrs[i]-1)>= ptrs[i+1]-1-nnzbase[count] 00123 if(mycoords[COL] == count){// ray belongs to count and it's a new starting ray 00124 nraysloc++; 00125 } 00126 ptrstart[count] = i; 00127 count++; 00128 } 00129 00130 } 00131 else { //ptrs[i]-1 > nnzbase[count] so ray belongs to count-1 00132 if(mycoords[COL] == count-1){ 00133 nraysloc++; 00134 } 00135 00136 } 00137 } // end for 00138 ptrstart[dims[COL]] = nrays; 00139 return nraysloc; 00140 00141 }
|
1.3.9.1