#include "mpi.h"
#include "stdlib.h"
#include "emdata.h"
#include "ali3d_unified_mpi.h"
#include "utilcomm.h"
Include dependency graph for utilcomm.cpp:
Go to the source code of this file.
Functions | |
int | ReadVandBcast (MPI_Comm comm, EMData *volume, char *volfname) |
int | ReadStackandDist (MPI_Comm comm, EMData ***images2D, char *stackfname, int *nloc) |
int | ReadAngTrandDist (MPI_Comm comm, float *angleshift, char *paramfname, int nloc) |
int | CleanStack (MPI_Comm comm, EMData **image_stack, int nloc, int ri, Vec3i volsize, Vec3i origin) |
int | setpart (MPI_Comm comm, int nang, int *psize, int *nbase) |
int CleanStack | ( | MPI_Comm | comm, | |
EMData ** | image_stack, | |||
int | nloc, | |||
int | ri, | |||
Vec3i | volsize, | |||
Vec3i | origin | |||
) |
Definition at line 233 of file utilcomm.cpp.
References asta2(), EMDeleteArray(), get_data(), imgdata, nx, ny, and rhs.
Referenced by main().
00234 { 00235 int nx = volsize[0]; 00236 int ny = volsize[1]; 00237 00238 float * rhs = new float[nx*ny]; 00239 float * imgdata; 00240 00241 // Calculate average "background" from all pixels strictly outside of radius 00242 double aba, abaloc; // average background 00243 aba = 0.0; 00244 abaloc = 0.0; 00245 int klp, klploc; // number of pixels in background 00246 klp = 0; 00247 klploc = 0; 00248 00249 // calculate avg background in parallel 00250 for ( int i = 0 ; i < nloc ; ++i ) { 00251 imgdata = image_stack[i]->get_data(); 00252 asta2(imgdata, nx, ny, ri, &abaloc, &klploc); 00253 } 00254 00255 MPI_Allreduce(&abaloc, &aba, 1, MPI_DOUBLE, MPI_SUM, comm); 00256 MPI_Allreduce(&klploc, &klp, 1, MPI_INT, MPI_SUM, comm); 00257 00258 aba /= klp; 00259 00260 // subtract off the average background from pixels weakly inside of radius 00261 int x_summand, y_summand; 00262 int r_squared = ri * ri; 00263 for ( int i = 0 ; i < nloc ; ++i ) { 00264 imgdata = image_stack[i]->get_data(); 00265 for ( int j = 0 ; j < nx ; ++j) { 00266 x_summand = (j - origin[0]) * (j - origin[0]); 00267 for ( int k = 0 ; k < ny ; ++k ) { 00268 y_summand = (k - origin[1]) * (k - origin[1]); 00269 if ( x_summand + y_summand <= r_squared) { 00270 imgdata[j*ny + k] -= aba; 00271 } 00272 } 00273 } 00274 } 00275 EMDeleteArray(rhs); 00276 return 0; 00277 }
int ReadAngTrandDist | ( | MPI_Comm | comm, | |
float * | angleshift, | |||
char * | paramfname, | |||
int | nloc | |||
) |
Definition at line 157 of file utilcomm.cpp.
References EMDeleteArray(), and ierr.
Referenced by main().
00158 { 00159 // read a parameter file that contains a list of angles and shifts and distribute 00160 // them evenly among different processors 00161 00162 int mypid, ierr=0, mpierr=0, nimgs=0, ncpus; 00163 double t0; 00164 MPI_Status mpistatus; 00165 FILE *fp=NULL; 00166 00167 float *iobuffer = new float[5*nloc]; 00168 if (!iobuffer) { 00169 fprintf(stderr,"failed to allocate buffer to read angles shifts\n"); 00170 ierr = -1; 00171 goto EXIT; 00172 } 00173 00174 MPI_Comm_rank(comm,&mypid); 00175 MPI_Comm_size(comm,&ncpus); 00176 00177 t0 = MPI_Wtime(); 00178 00179 if (mypid ==0) { 00180 fp = fopen(paramfname,"r"); 00181 if (!fp) ierr = 1; 00182 } 00183 MPI_Bcast(&ierr, 1, MPI_INT, 0, comm); 00184 00185 if ( ierr ) { 00186 if (mypid ==0) fprintf(stderr,"failed to open %s\n", paramfname); 00187 } 00188 else { 00189 if (mypid == 0) { 00190 for (int iproc = 0; iproc < ncpus; iproc++) { 00191 // figure out the number of images assigned to processor iproc 00192 if (iproc > 0) { 00193 MPI_Recv(&nimgs, 1, MPI_INT, iproc, iproc, comm, &mpistatus); 00194 // Read the next nimgs set of angles and shifts 00195 for (int i = 0; i < nimgs; i++) { 00196 fscanf(fp,"%f %f %f %f %f", 00197 &iobuffer[5*i+0], 00198 &iobuffer[5*i+1], 00199 &iobuffer[5*i+2], 00200 &iobuffer[5*i+3], 00201 &iobuffer[5*i+4]); 00202 } 00203 MPI_Send(iobuffer,5*nimgs,MPI_FLOAT,iproc,iproc,comm); 00204 } 00205 else { 00206 for (int i = 0; i < nloc; i++) { 00207 fscanf(fp,"%f %f %f %f %f", 00208 &angleshift[5*i+0], 00209 &angleshift[5*i+1], 00210 &angleshift[5*i+2], 00211 &angleshift[5*i+3], 00212 &angleshift[5*i+4]); 00213 } 00214 } 00215 } 00216 fclose(fp); 00217 } 00218 else { 00219 // send image count to the master processor (mypid = 0) 00220 MPI_Send(&nloc, 1, MPI_INT, 0, mypid, comm); 00221 // Receive angleshifts 00222 MPI_Recv(angleshift, 5*nloc, MPI_FLOAT, 0, mypid, comm, &mpistatus); 00223 } 00224 } 00225 EMDeleteArray(iobuffer); 00226 if (mypid == 0) 00227 printf("I/O time for reading angles & shifts = %11.3e\n", 00228 MPI_Wtime() - t0); 00229 EXIT: 00230 return ierr; 00231 }
int ReadStackandDist | ( | MPI_Comm | comm, | |
EMData *** | images2D, | |||
char * | stackfname, | |||
int * | nloc | |||
) |
Definition at line 48 of file utilcomm.cpp.
References EMDeleteArray(), EMDeletePtr(), EMAN::EMUtil::get_image_count(), ierr, img_ptr, imgdata, nx, ny, and setpart().
Referenced by main().
00049 { 00050 int ncpus, mypid, ierr, mpierr=0; 00051 MPI_Status mpistatus; 00052 EMUtil *my_util; 00053 int nima; 00054 FILE *fp=NULL; 00055 00056 MPI_Comm_size(comm,&ncpus); 00057 MPI_Comm_rank(comm,&mypid); 00058 00059 ierr = 0; 00060 if (mypid == 0) { 00061 fp = fopen(stackfname,"r"); 00062 if (!fp) { 00063 ierr = 1; 00064 printf("failed to open %s\n", stackfname); 00065 } 00066 else { 00067 fclose(fp); 00068 nima = my_util->get_image_count(stackfname); 00069 } 00070 } 00071 mpierr = MPI_Bcast(&ierr, 1, MPI_INT, 0, comm); 00072 if (ierr == 0) { 00073 int *psize = new int[ncpus]; 00074 int *nbase = new int[ncpus]; 00075 00076 mpierr = MPI_Bcast(&nima, 1, MPI_INT, 0, comm); 00077 00078 *nloc = setpart(comm, nima, psize, nbase); 00079 *images2D = new EMData*[*nloc]; // NB!: whoever calls ReadStackandDist must delete this! 00080 00081 EMData *img_ptr; 00082 int img_index; 00083 float *imgdata; 00084 00085 // read the first image to get size 00086 img_ptr = new EMData(); 00087 img_ptr->read_image(stackfname, 0); 00088 int nx = img_ptr->get_xsize(); 00089 int ny = img_ptr->get_ysize(); 00090 00091 float s2x, s2y; 00092 00093 if (mypid == 0) { 00094 printf("Master node reading and distributing %d images...\n", nima); 00095 for ( int ip = 0 ; ip < ncpus ; ++ip ) { 00096 for ( int i = 0 ; i < psize[ip] ; ++i ) { 00097 img_index = nbase[ip] + i; 00098 if (ip != 0) { 00099 img_ptr->read_image(stackfname, img_index); 00100 // get a pointer to the image's data 00101 imgdata = img_ptr->get_data(); 00102 // find the x/y shift values if it has them, otherwise set them to 0.0 00103 try { 00104 s2x = (*images2D)[i]->get_attr("s2x"); 00105 } catch ( std::exception& e ) { 00106 s2x = 0.0; 00107 } 00108 try { 00109 s2y = (*images2D)[i]->get_attr("s2y"); 00110 } catch ( std::exception& e ) { 00111 s2y = 0.0; 00112 } 00113 // send these to processor ip 00114 MPI_Send(imgdata, nx*ny, MPI_FLOAT, ip, ip, comm); 00115 MPI_Send(&s2x, 1, MPI_FLOAT, ip, ip, comm); 00116 MPI_Send(&s2y, 1, MPI_FLOAT, ip, ip, comm); 00117 } else { // ip == 0 00118 (*images2D)[i] = new EMData(); 00119 (*images2D)[i]->read_image(stackfname, img_index); 00120 try { 00121 s2x = (*images2D)[i]->get_attr("s2x"); 00122 } catch ( std::exception& e ) { 00123 (*images2D)[i]->set_attr("s2x",0.0); 00124 } 00125 try { 00126 s2y = (*images2D)[i]->get_attr("s2y"); 00127 } catch ( std::exception& e ) { 00128 (*images2D)[i]->set_attr("s2y",0.0); 00129 } 00130 } 00131 } 00132 printf("finished reading data for processor %d\n", ip); 00133 } 00134 } else { // mypid != 0 : everyone else receives and reads in their data 00135 for ( int i = 0 ; i < psize[mypid] ; ++i ) { 00136 (*images2D)[i] = new EMData(); 00137 (*images2D)[i]->set_size(nx, ny, 1); 00138 imgdata = (*images2D)[i]->get_data(); 00139 MPI_Recv(imgdata, nx*ny, MPI_FLOAT, 0, mypid, comm, &mpistatus); 00140 MPI_Recv(&s2x, 1, MPI_FLOAT, 0, mypid, comm, &mpistatus); 00141 MPI_Recv(&s2y, 1, MPI_FLOAT, 0, mypid, comm, &mpistatus); 00142 (*images2D)[i]->set_attr("s2x",s2x); 00143 (*images2D)[i]->set_attr("s2y",s2y); 00144 } 00145 printf("received %d images for processor %d\n", *nloc, mypid); 00146 } 00147 if (mypid == 0) printf("finished reading and distributing data\n"); 00148 ierr = mpierr; 00149 00150 EMDeletePtr(img_ptr); 00151 EMDeleteArray(psize); 00152 EMDeleteArray(nbase); 00153 } 00154 return ierr; 00155 }
int ReadVandBcast | ( | MPI_Comm | comm, | |
EMData * | volume, | |||
char * | volfname | |||
) |
Definition at line 10 of file utilcomm.cpp.
References EMAN::EMData::get_data(), EMAN::EMData::get_ndim(), EMAN::EMData::get_xsize(), EMAN::EMData::get_ysize(), EMAN::EMData::get_zsize(), ierr, nx, ny, EMAN::EMData::read_image(), and EMAN::EMData::set_size().
Referenced by main().
00011 { 00012 int mypid, ierr, ndim, nx, ny, nz, mpierr=0; 00013 MPI_Comm_rank(comm,&mypid); 00014 FILE *fp=NULL; 00015 00016 ierr = 0; 00017 if (mypid == 0) { 00018 // check to see if the file exists 00019 fp = fopen(volfname,"r"); 00020 if (!fp) { 00021 ierr = 1; 00022 printf("failed to open %s\n", volfname); 00023 } 00024 else { 00025 fclose(fp); 00026 volume->read_image(volfname); 00027 ndim = volume->get_ndim(); 00028 nx = volume->get_xsize(); 00029 ny = volume->get_ysize(); 00030 nz = volume->get_zsize(); 00031 printf("ndim = %d, nx = %d, ny = %d, nz = %d\n", ndim, nx, ny, nz); 00032 } 00033 } 00034 mpierr = MPI_Bcast(&ierr, 1, MPI_INT, 0, comm); 00035 if (ierr==0) { 00036 mpierr = MPI_Bcast(&nx, 1, MPI_INT, 0, comm); 00037 mpierr = MPI_Bcast(&ny, 1, MPI_INT, 0, comm); 00038 mpierr = MPI_Bcast(&nz, 1, MPI_INT, 0, comm); 00039 if (mypid !=0) volume->set_size(nx,ny,nz); 00040 00041 float * voldata = volume->get_data(); 00042 mpierr = MPI_Bcast(voldata, nx*ny*nz, MPI_FLOAT, 0, comm); 00043 ierr = mpierr; 00044 } 00045 return ierr; 00046 }
int setpart | ( | MPI_Comm | comm, | |
int | nang, | |||
int * | psize, | |||
int * | nbase | |||
) |
Definition at line 279 of file utilcomm.cpp.
Referenced by ali3d_d(), fgcalc(), ReadStackandDist(), and unified().
00280 { 00281 int ncpus, mypid, nangloc, nrem; 00282 00283 MPI_Comm_size(comm,&ncpus); 00284 MPI_Comm_rank(comm,&mypid); 00285 00286 nangloc = nang/ncpus; 00287 nrem = nang - ncpus*nangloc; 00288 if (mypid < nrem) nangloc++; 00289 00290 for (int i = 0; i < ncpus; i++) { 00291 psize[i] = nang/ncpus; 00292 if (i < nrem) psize[i]++; 00293 } 00294 00295 nbase[0] = 0; 00296 for (int i = 1; i < ncpus; i++) { 00297 nbase[i] = nbase[i-1] + psize[i-1]; 00298 } 00299 00300 return nangloc; 00301 }