#include <stdio.h>
#include <stdlib.h>
#include <mpi.h>
#include <math.h>
#include <string.h>
Include dependency graph for runcartrec.cpp:
Go to the source code of this file.
Functions | |
int | setpart_gc1 (MPI_Comm comm_2d, int nangs, int *psize, int *nbase) |
int | getnnz (int *volsize, int ri, int *origin, int *nrays, int *nnz) |
int | setpart_gr1 (MPI_Comm comm_2d, int nnz, int *nnzpart, int *nnzbase) |
int | getcb2sph (int *volsize, int ri, int *origin, int nnz0, int *ptrs, int *cord) |
int | sphpart (MPI_Comm comm_2d, int nrays, int *ptrs, int *nnzbase, int *ptrstart) |
int | make_proj_mat (float phi, float theta, float psi, float *dm) |
int | ifix (float a) |
int | fwdpj3_Cart (int *volsize, int nraysloc, int nnzloc, float *dm, int *origin, int ri, int *ptrs, int *cord, int myptrstart, float *x, float *y) |
int | bckpj3_Cart (int *volsize, int nraysloc, int nnzloc, float *dm, int *origin, int ri, int *ptrs, int *cord, int myptrstart, float *x, float *y) |
int | fwdpj3 (int *volsize, int nrays, int nnz, float *dm, int *origin, int ri, int *ptrs, int *cord, float *x, float *y) |
int | bckpj3 (int *volsize, int nrays, int nnz, float *dm, int *origin, int ri, int *ptrs, int *cord, float *x, float *y) |
int | main (int argc, char *argv[]) |
Variables | |
int | nangs |
int | nangsloc |
int | nx |
int | ny |
int | myptrstart |
int | nnz |
int | nrays |
int | nnzloc |
int | nraysloc |
float * | imagesloc |
float * | images |
float * | newimagesloc |
float * | reprojloc |
float * | angles |
float * | anglesloc |
float * | vol_sph |
float * | onevol_sph |
float * | reproj |
int bckpj3 | ( | int * | volsize, | |
int | nrays, | |||
int | nnz, | |||
float * | dm, | |||
int * | origin, | |||
int | ri, | |||
int * | ptrs, | |||
int * | cord, | |||
float * | x, | |||
float * | y | |||
) |
Definition at line 618 of file utilcomm2d.cpp.
References cord, dm, ifix(), nx, ny, ptrs, status, x, and y.
00619 { 00620 int i, j, iqx,iqy, xc, yc, zc; 00621 float xb, yb, dx, dy, dx1m, dy1m, dxdy; 00622 int status = 0; 00623 00624 int xcent = origin[0]; 00625 int ycent = origin[1]; 00626 int zcent = origin[2]; 00627 00628 int nx = volsize[0]; 00629 int ny = volsize[1]; 00630 00631 // Phi: adding the shift parameters that get passed in as the last two entries of dm 00632 float sx, sy; 00633 00634 sx = dm(7); 00635 sy = dm(8); 00636 00637 00638 if ( nx > 2*ri) { 00639 for (i = 1; i <= nrays; i++) { 00640 zc = cord(1,i) - zcent; 00641 yc = cord(2,i) - ycent; 00642 xc = cord(3,i) - xcent; 00643 00644 xb = zc*dm(1)+yc*dm(2)+xc*dm(3) + xcent + sx; 00645 yb = zc*dm(4)+yc*dm(5)+xc*dm(6) + ycent + sy; 00646 00647 for (j = ptrs(i); j <ptrs(i+1); j++) { 00648 iqx = ifix(xb); 00649 iqy = ifix(yb); 00650 00651 dx = xb - iqx; 00652 dy = yb - iqy; 00653 dx1m = 1.0 - dx; 00654 dy1m = 1.0 - dy; 00655 dxdy = dx*dy; 00656 /* 00657 c y(j) = y(j) + dx1m*dy1m*x(iqx , iqy) 00658 c & + dx1m*dy *x(iqx , iqy+1) 00659 c & + dx *dy1m*x(iqx+1, iqy) 00660 c & + dx *dy *x(iqx+1, iqy+1) 00661 c 00662 c --- faster version of the above commented out 00663 c code (derived by summing the following table 00664 c of coefficients along the colunms) --- 00665 c 00666 c 1 dx dy dxdy 00667 c ------ -------- -------- ------- 00668 c x(i,j) -x(i,j) -x(i,j) x(i,j) 00669 c x(i,j+1) -x(i,j+1) 00670 c x(i+1,j) -x(i+1,j) 00671 c x(i+1,j+1) 00672 c 00673 */ 00674 // Phi: add index checking, now that shifts are being used 00675 if ( iqx <= nx && iqy <= ny && iqx >= 1 && iqy >= 1 ) { 00676 y(j) += x(iqx,iqy); 00677 if ( iqx + 1 <= nx && iqx + 1 >= 1 ) { 00678 y(j) += dx*(-x(iqx,iqy)+x(iqx+1,iqy)); 00679 } 00680 if ( iqy + 1 <= ny && iqy + 1 >= 1 ) { 00681 y(j) += dy*(-x(iqx,iqy)+x(iqx,iqy+1)); 00682 } 00683 if ( iqx + 1 <= nx && iqy + 1 <= ny && iqx + 1 >= 1 && iqy + 1 >= 1 ) { 00684 y(j) += dxdy*( x(iqx,iqy) - x(iqx,iqy+1) -x(iqx+1,iqy) + x(iqx+1,iqy+1) ); 00685 } 00686 } 00687 00688 // y(j) += x(iqx,iqy) 00689 // + dx*(-x(iqx,iqy)+x(iqx+1,iqy)) 00690 // + dy*(-x(iqx,iqy)+x(iqx,iqy+1)) 00691 // + dxdy*( x(iqx,iqy) - x(iqx,iqy+1) 00692 // -x(iqx+1,iqy) + x(iqx+1,iqy+1) ); 00693 00694 xb += dm(1); 00695 yb += dm(4); 00696 } // end for j 00697 } // end for i 00698 } 00699 else { 00700 fprintf(stderr, "bckpj3: nx must be greater than 2*ri\n"); 00701 } 00702 00703 return status; 00704 }
int bckpj3_Cart | ( | int * | volsize, | |
int | nraysloc, | |||
int | nnzloc, | |||
float * | dm, | |||
int * | origin, | |||
int | ri, | |||
int * | ptrs, | |||
int * | cord, | |||
int | myptrstart, | |||
float * | x, | |||
float * | y | |||
) |
Definition at line 412 of file utilcomm2d.cpp.
References cord, dm, ifix(), nx, ny, ptrs, status, x, and y.
00415 : 00416 volsize - size of 3D cube volume 00417 nraysloc - local number of rays 00418 nnzloc - local number of voxels 00419 dm - projection matrix 00420 origin, ri - origin and radius of sphere 00421 ptrs, cord - pointers and coordinates for first ray 00422 myptrstart - determines starting rays 00423 x - 2D image 00424 Output: 00425 y - portion of the 3D volume 00426 */ 00427 { 00428 int i, j, iqx,iqy, xc, yc, zc, jj; 00429 float xb, yb, dx, dy, dx1m, dy1m, dxdy; 00430 int status = 0; 00431 00432 int xcent = origin[0]; 00433 int ycent = origin[1]; 00434 int zcent = origin[2]; 00435 00436 int nx = volsize[0]; 00437 int ny = volsize[1]; 00438 00439 // Phi: adding the shift parameters that get passed in as the last two entries of dm 00440 float sx, sy; 00441 00442 sx = dm(7); 00443 sy = dm(8); 00444 00445 if ( nx > 2*ri) { 00446 for (i = 1; i <= nraysloc; i++) { 00447 zc = cord(1,myptrstart+i) - zcent; 00448 yc = cord(2,myptrstart+i) - ycent; 00449 xc = cord(3,myptrstart+i) - xcent; 00450 00451 xb = zc*dm(1)+yc*dm(2)+xc*dm(3) + xcent + sx; 00452 yb = zc*dm(4)+yc*dm(5)+xc*dm(6) + ycent + sy; 00453 00454 for (j = ptrs(myptrstart+i); j <ptrs(myptrstart+i+1); j++) { 00455 jj = j-ptrs[myptrstart]+1; 00456 00457 // jj is the index on the local volume 00458 iqx = ifix(xb); 00459 iqy = ifix(yb); 00460 00461 dx = xb - iqx; 00462 dy = yb - iqy; 00463 dx1m = 1.0 - dx; 00464 dy1m = 1.0 - dy; 00465 dxdy = dx*dy; 00466 /* 00467 c y(j) = y(j) + dx1m*dy1m*x(iqx , iqy) 00468 c & + dx1m*dy *x(iqx , iqy+1) 00469 c & + dx *dy1m*x(iqx+1, iqy) 00470 c & + dx *dy *x(iqx+1, iqy+1) 00471 c 00472 c --- faster version of the above commented out 00473 c code (derived by summing the following table 00474 c of coefficients along the colunms) --- 00475 c 00476 c 1 dx dy dxdy 00477 c ------ -------- -------- ------- 00478 c x(i,j) -x(i,j) -x(i,j) x(i,j) 00479 c x(i,j+1) -x(i,j+1) 00480 c x(i+1,j) -x(i+1,j) 00481 c x(i+1,j+1) 00482 c 00483 */ 00484 // Phi: add index checking, now that shifts are being used 00485 if ( iqx <= nx && iqy <= ny && iqx >= 1 && iqy >= 1 ) { 00486 y(jj) += x(iqx,iqy); 00487 if ( iqx + 1 <= nx && iqx + 1 >= 1 ) { 00488 y(jj) += dx*(-x(iqx,iqy)+x(iqx+1,iqy)); 00489 } 00490 if ( iqy + 1 <= ny && iqy + 1 >= 1 ) { 00491 y(jj) += dy*(-x(iqx,iqy)+x(iqx,iqy+1)); 00492 } 00493 if ( iqx + 1 <= nx && iqy + 1 <= ny && iqx + 1 >= 1 && iqy + 1 >= 1 ) { 00494 y(jj) += dxdy*( x(iqx,iqy) - x(iqx,iqy+1) -x(iqx+1,iqy) + x(iqx+1,iqy+1) ); 00495 } 00496 } 00497 00498 // y(j) += x(iqx,iqy) 00499 // + dx*(-x(iqx,iqy)+x(iqx+1,iqy)) 00500 // + dy*(-x(iqx,iqy)+x(iqx,iqy+1)) 00501 // + dxdy*( x(iqx,iqy) - x(iqx,iqy+1) 00502 // -x(iqx+1,iqy) + x(iqx+1,iqy+1) ); 00503 00504 xb += dm(1); 00505 yb += dm(4); 00506 } // end for j 00507 } // end for i 00508 } 00509 else { 00510 fprintf(stderr, "bckpj3: nx must be greater than 2*ri\n"); 00511 } 00512 00513 return status; 00514 }
int fwdpj3 | ( | int * | volsize, | |
int | nrays, | |||
int | nnz, | |||
float * | dm, | |||
int * | origin, | |||
int | ri, | |||
int * | ptrs, | |||
int * | cord, | |||
float * | x, | |||
float * | y | |||
) |
Definition at line 526 of file utilcomm2d.cpp.
References cord, dm, ifix(), nx, ny, ptrs, status, x, and y.
00527 { 00528 /* 00529 purpose: y <--- proj(x) 00530 input : volsize the size (nx,ny,nz) of the volume 00531 nrays number of rays within the compact spherical 00532 representation 00533 nnz number of voxels within the sphere 00534 dm an array of size 9 storing transformation 00535 associated with the projection direction 00536 origin coordinates of the center of the volume 00537 ri radius of the sphere 00538 ptrs the beginning address of each ray 00539 cord the coordinates of the first point in each ray 00540 x 3d input volume 00541 y 2d output image 00542 */ 00543 00544 int iqx, iqy, i, j, xc, yc, zc; 00545 float ct, dipx, dipy, dipx1m, dipy1m, xb, yb, dm1, dm4; 00546 int status = 0; 00547 00548 // Phi: adding the shift parameters that get passed in as the last two entries of dm 00549 float sx, sy; 00550 00551 sx = dm(7); 00552 sy = dm(8); 00553 00554 int xcent = origin[0]; 00555 int ycent = origin[1]; 00556 int zcent = origin[2]; 00557 00558 int nx = volsize[0]; 00559 int ny = volsize[1]; 00560 00561 dm1 = dm(1); 00562 dm4 = dm(4); 00563 00564 if ( nx > 2*ri ) { 00565 for (i = 1; i <= nrays; i++) { 00566 00567 zc = cord(1,i)-zcent; 00568 yc = cord(2,i)-ycent; 00569 xc = cord(3,i)-xcent; 00570 xb = zc* dm(1) +yc* dm(2) +xc* dm(3) + xcent + sx; 00571 yb = zc* dm(4) +yc* dm(5) +xc* dm(6) + ycent + sy; 00572 00573 for (j = ptrs(i); j< ptrs(i+1); j++) { 00574 iqx = ifix(xb); 00575 iqy = ifix(yb); 00576 00577 ct = x(j); 00578 00579 // dipx = xb - (float)(iqx); 00580 // dipy = (yb - (float)(iqy)) * ct; 00581 dipx = xb - iqx; 00582 dipy = (yb - iqy) * ct; 00583 00584 dipy1m = ct - dipy; 00585 dipx1m = 1.0 - dipx; 00586 00587 if (iqx <= nx && iqy <= ny && iqx >= 1 && iqy >= 1) 00588 // y(iqx ,iqy) = y(iqx ,iqy) + dipx1m*dipy1m; 00589 y(iqx ,iqy) += dipx1m*dipy1m; 00590 if (iqx + 1 <= nx && iqy <= ny && iqx >= 0 && iqy >= 1) 00591 // y(iqx+1,iqy) = y(iqx+1,iqy) + dipx*dipy1m; 00592 y(iqx+1,iqy) += dipx*dipy1m; 00593 if (iqx + 1 <= nx && iqy + 1 <= ny && iqx >= 0 && iqy >= 0) 00594 // y(iqx+1,iqy+1) = y(iqx+1,iqy+1) + dipx*dipy; 00595 y(iqx+1,iqy+1) += dipx*dipy; 00596 if (iqx <= nx && iqy + 1 <= ny && iqx >= 1 && iqy >= 0) 00597 // y(iqx ,iqy+1) = y(iqx ,iqy+1) + dipx1m*dipy; 00598 y(iqx ,iqy+1) += dipx1m*dipy; 00599 xb += dm1; 00600 yb += dm4; 00601 } 00602 } 00603 } 00604 else { 00605 fprintf(stderr, " nx must be greater than 2*ri\n"); 00606 exit(1); 00607 } 00608 return status; 00609 }
int fwdpj3_Cart | ( | int * | volsize, | |
int | nraysloc, | |||
int | nnzloc, | |||
float * | dm, | |||
int * | origin, | |||
int | ri, | |||
int * | ptrs, | |||
int * | cord, | |||
int | myptrstart, | |||
float * | x, | |||
float * | y | |||
) |
Definition at line 322 of file utilcomm2d.cpp.
References cord, dm, ifix(), nx, ny, ptrs, status, x, and y.
00325 : 00326 volsize - size (nx,ny,nz) of 3D cube volume 00327 nraysloc - local number of rays within the compact spherical representation 00328 nnzloc - local number of voxels 00329 dm - an array of size 9 storing transformation 00330 origin, ri - origin and radius of sphere 00331 ptrs, cord - pointers and coordinates for first ray 00332 myptrstart - determines starting rays 00333 x - portion of the 3D input volume 00334 Output: 00335 y - 2D output image 00336 */ 00337 { 00338 int iqx, iqy, i, j, xc, yc, zc, jj; 00339 float ct, dipx, dipy, dipx1m, dipy1m, xb, yb, dm1, dm4; 00340 int status = 0; 00341 00342 // Phi: adding the shift parameters that get passed in as the last two entries of dm 00343 float sx, sy; 00344 00345 sx = dm(7); 00346 sy = dm(8); 00347 00348 int xcent = origin[0]; 00349 int ycent = origin[1]; 00350 int zcent = origin[2]; 00351 00352 int nx = volsize[0]; 00353 int ny = volsize[1]; 00354 00355 dm1 = dm(1); 00356 dm4 = dm(4); 00357 00358 if ( nx > 2*ri ) { 00359 for (i = 1; i <= nraysloc; i++) { 00360 00361 zc = cord(1,myptrstart+i)-zcent; 00362 yc = cord(2,myptrstart+i)-ycent; 00363 xc = cord(3,myptrstart+i)-xcent; 00364 xb = zc* dm(1) +yc* dm(2) +xc* dm(3) + xcent + sx; 00365 yb = zc* dm(4) +yc* dm(5) +xc* dm(6) + ycent + sy; 00366 00367 for (j = ptrs(myptrstart+i); j< ptrs(myptrstart+i+1); j++) { 00368 jj = j-ptrs[myptrstart]+1; 00369 iqx = ifix(xb); 00370 iqy = ifix(yb); 00371 00372 // jj is the index on the local volume 00373 ct = x(jj); 00374 00375 // dipx = xb - (float)(iqx); 00376 // dipy = (yb - (float)(iqy)) * ct; 00377 dipx = xb - iqx; 00378 dipy = (yb - iqy) * ct; 00379 00380 dipy1m = ct - dipy; 00381 dipx1m = 1.0 - dipx; 00382 00383 if (iqx <= nx && iqy <= ny && iqx >= 1 && iqy >= 1) 00384 // y(iqx ,iqy) = y(iqx ,iqy) + dipx1m*dipy1m; 00385 y(iqx ,iqy) += dipx1m*dipy1m; 00386 if (iqx + 1 <= nx && iqy <= ny && iqx >= 0 && iqy >= 1) 00387 // y(iqx+1,iqy) = y(iqx+1,iqy) + dipx*dipy1m; 00388 y(iqx+1,iqy) += dipx*dipy1m; 00389 if (iqx + 1 <= nx && iqy + 1 <= ny && iqx >= 0 && iqy >= 0) 00390 // y(iqx+1,iqy+1) = y(iqx+1,iqy+1) + dipx*dipy; 00391 y(iqx+1,iqy+1) += dipx*dipy; 00392 if (iqx <= nx && iqy + 1 <= ny && iqx >= 1 && iqy >= 0) 00393 // y(iqx ,iqy+1) = y(iqx ,iqy+1) + dipx1m*dipy; 00394 y(iqx ,iqy+1) += dipx1m*dipy; 00395 xb += dm1; 00396 yb += dm4; 00397 } 00398 } 00399 } 00400 else { 00401 fprintf(stderr, " nx must be greater than 2*ri\n"); 00402 exit(1); 00403 } 00404 return status; 00405 } #undef x
int getcb2sph | ( | int * | volsize, | |
int | ri, | |||
int * | origin, | |||
int | nnz0, | |||
int * | ptrs, | |||
int * | cord | |||
) |
Definition at line 144 of file utilcomm2d.cpp.
References cord, nnz, nrays, nx, ny, ptrs, and status.
00147 : 00148 volsize - dimensions of 3D cube 00149 ri - radius of sphere 00150 origin - origin of sphere 00151 nnz0 - number of non-zeros 00152 00153 Output: 00154 ptrs - vector containing pointers to the beginning of each ray 00155 cord - vector containing the x,y,z coordinates of the beginning of each ray 00156 */ 00157 00158 { 00159 int xs, ys, zs, xx, yy, zz, rs, r2; 00160 int ix, iy, iz, jnz, nnz, nrays; 00161 int ftm = 0, status = 0; 00162 00163 int xcent = (int)origin[0]; 00164 int ycent = (int)origin[1]; 00165 int zcent = (int)origin[2]; 00166 00167 int nx = (int)volsize[0]; 00168 int ny = (int)volsize[1]; 00169 int nz = (int)volsize[2]; 00170 00171 r2 = ri*ri; 00172 nnz = 0; 00173 nrays = 0; 00174 ptrs(1) = 1; 00175 00176 for (ix = 1; ix <= nx; ix++) { 00177 xs = ix-xcent; 00178 xx = xs*xs; 00179 for ( iy = 1; iy <= ny; iy++ ) { 00180 ys = iy-ycent; 00181 yy = ys*ys; 00182 jnz = 0; 00183 00184 ftm = 1; 00185 // not the most efficient implementation 00186 for (iz = 1; iz <= nz; iz++) { 00187 zs = iz-zcent; 00188 zz = zs*zs; 00189 rs = xx + yy + zz; 00190 if (rs <= r2) { 00191 jnz++; 00192 nnz++; 00193 // sphere(nnz) = cube(iz, iy, ix); 00194 00195 // record the coordinates of the first nonzero === 00196 if (ftm) { 00197 nrays++; 00198 cord(1,nrays) = iz; 00199 cord(2,nrays) = iy; 00200 cord(3,nrays) = ix; 00201 ftm = 0; 00202 } 00203 } 00204 } // end for (iz..) 00205 if (jnz > 0) { 00206 ptrs(nrays+1) = ptrs(nrays) + jnz; 00207 } // endif (jnz) 00208 } // end for iy 00209 } // end for ix 00210 if (nnz != nnz0) status = -1; 00211 return status; 00212 }
int getnnz | ( | int * | volsize, | |
int | ri, | |||
int * | origin, | |||
int * | nrays, | |||
int * | nnz | |||
) |
Definition at line 43 of file utilcomm2d.cpp.
References nx, ny, and status.
00045 : count the number of voxels within a sphere centered 00046 at origin and with a radius ri. 00047 00048 input: 00049 volsize vector containing size of the volume/cube 00050 ri radius of the object embedded in the cube. 00051 origin coordinates for the center of the volume 00052 00053 output: 00054 nnz total number of voxels within the sphere (of radius ri) 00055 nrays number of rays in z-direction. 00056 */ 00057 { 00058 int ix, iy, iz, rs, r2, xs, ys, zs, xx, yy, zz; 00059 int ftm=0, status = 0; 00060 00061 r2 = ri*ri; 00062 *nnz = 0; 00063 *nrays = 0; 00064 int nx = volsize[0]; 00065 int ny = volsize[1]; 00066 int nz = volsize[2]; 00067 // int nx = (int)volsize[0]; 00068 // int ny = (int)volsize[1]; 00069 // int nz = (int)volsize[2]; 00070 00071 int xcent = origin[0]; 00072 int ycent = origin[1]; 00073 int zcent = origin[2]; 00074 // int xcent = (int)origin[0]; 00075 // int ycent = (int)origin[1]; 00076 // int zcent = (int)origin[2]; 00077 00078 // need to add some error checking 00079 for (ix = 1; ix <=nx; ix++) { 00080 xs = ix-xcent; 00081 xx = xs*xs; 00082 for (iy = 1; iy <= ny; iy++) { 00083 ys = iy-ycent; 00084 yy = ys*ys; 00085 ftm = 1; 00086 for (iz = 1; iz <= nz; iz++) { 00087 zs = iz-zcent; 00088 zz = zs*zs; 00089 rs = xx + yy + zz; 00090 if (rs <= r2) { 00091 (*nnz)++; 00092 if (ftm) { 00093 (*nrays)++; 00094 ftm = 0; 00095 } 00096 } 00097 } 00098 } // end for iy 00099 } // end for ix 00100 return status; 00101 }
int ifix | ( | float | a | ) |
Definition at line 41 of file utilnum.cpp.
00042 { 00043 int ia = 0; 00044 if (a >= 0.0) { 00045 ia = (int)floor(a); 00046 } 00047 else { 00048 ia = (int)ceil(a); 00049 } 00050 return ia; 00051 }
int main | ( | int | argc, | |
char * | argv[] | |||
) |
Definition at line 32 of file runcartrec.cpp.
References angles, anglesloc, bckpj3(), bckpj3_Cart(), COL, cord, dm, fwdpj3(), fwdpj3_Cart(), getcb2sph(), getnnz(), ierr, images, imagesloc, make_proj_mat(), max, myptrstart, nangs, nangsloc, newimagesloc, nnz, nnzloc, nrays, nraysloc, nx, ny, onevol_sph, phi, ptrs, reproj, reprojloc, ROW, setpart_gc1(), setpart_gr1(), sphpart(), sqrt(), theta, and vol_sph.
00033 { 00034 int ncpus, mypid, nrem, ierr; 00035 MPI_Status mpistatus; 00036 MPI_Comm comm = MPI_COMM_WORLD; 00037 00038 // Variables needed for the Cartesian topology. 00039 int ROW = 0, COL = 1; 00040 int dims[2], periods[2], keep_dims[2]; 00041 int my2dpid, mycoords[2], srcoords[2], otherpid; 00042 MPI_Comm comm_2d, comm_row, comm_col; 00043 00044 // Initialize MPI. 00045 MPI_Init(&argc, &argv); 00046 MPI_Comm_size(comm, &ncpus); 00047 MPI_Comm_rank(comm, &mypid); 00048 00049 if ( argc < 3 ) { 00050 printf ("ERROR: %s requires Cartesian dimensions input\n", argv[0]); 00051 return -1; 00052 } 00053 // Set up a Cartesian virtual topology and get the rank and coordinates of the processes in the topology. 00054 dims[ROW] = atoi(argv[1]); // Row dimension of the topology 00055 dims[COL] = atoi(argv[2]); // Column dimension of the topology 00056 00057 if (dims[ROW]*dims[COL] != ncpus){ 00058 printf("ERROR: Row dim and col dim not equal to ncpus\n"); 00059 return -1; 00060 } 00061 00062 periods[ROW] = periods[COL] = 1; // Set the periods for wrap-around 00063 00064 MPI_Cart_create(comm, 2, dims, periods, 1, &comm_2d); 00065 MPI_Comm_rank(comm_2d, &my2dpid); //Get my pid in the new 2D topology 00066 MPI_Cart_coords(comm_2d, my2dpid, 2, mycoords); // Get my coordinates 00067 00068 /* Create the row-based sub-topology */ 00069 keep_dims[ROW] = 0; 00070 keep_dims[COL] = 1; 00071 MPI_Cart_sub(comm_2d, keep_dims, &comm_row); 00072 00073 /* Create the column-based sub-topology */ 00074 keep_dims[ROW] = 1; 00075 keep_dims[COL] = 0; 00076 MPI_Cart_sub(comm_2d, keep_dims, &comm_col); 00077 00078 // STEP 1: Have processor (0,0) read in the entire set of 2D images, divide up the images, and send corresponding images to processors in processor group: g_c_0 Do the same for the angles. 00079 00080 if (mycoords[ROW] == 0 && mycoords[COL] == 0){ //I'm processor (0,0) 00081 FILE *fp, *fpa; 00082 char imagefname[80]="tf2d84.raw", anglesfname[80]="angles.dat"; 00083 00084 fp = fopen(imagefname,"r"); 00085 fread(&nangs, sizeof(int), 1, fp); 00086 fread(&nx, sizeof(int), 1, fp); 00087 fread(&ny, sizeof(int), 1, fp); 00088 00089 images = new float[nx*ny*nangs]; 00090 fread(images, sizeof(float), nx*ny*nangs, fp); 00091 fclose(fp); 00092 00093 fpa = fopen(anglesfname,"r"); 00094 angles = new float[3*nangs]; 00095 for (int i = 0; i< 3*nangs; i++) 00096 fscanf(fpa, "%f",&angles[i]); 00097 00098 fclose(fpa); 00099 printf("There are %d 2D images of size %d x %d\n", nangs, nx, ny); 00100 } 00101 00102 // Broadcast variables nangs, nx, ny to all processors 00103 srcoords[ROW] = srcoords[COL] = 0; 00104 MPI_Cart_rank(comm_2d, srcoords, &otherpid); 00105 00106 MPI_Bcast (&nangs, 1, MPI_INT, otherpid, comm_2d); 00107 MPI_Bcast (&nx, 1, MPI_INT, otherpid, comm_2d); 00108 MPI_Bcast (&ny, 1, MPI_INT, otherpid, comm_2d); 00109 00110 // Send images and angles from Processor (0,0) to processors in group g_c_0 00111 int *psize = new int[dims[ROW]]; 00112 int *nbase = new int[dims[ROW]]; 00113 00114 nangsloc = setpart_gc1(comm_2d, nangs, psize, nbase); 00115 imagesloc = new float[psize[mycoords[ROW]]*nx*ny]; 00116 reprojloc = new float[psize[mycoords[ROW]]*nx*ny]; 00117 anglesloc = new float[psize[mycoords[ROW]]*3]; 00118 00119 // printf("My coords are (%d,%d) and nangsloc = %d\n", mycoords[ROW], mycoords[COL], nangsloc); 00120 00121 if (mycoords[COL] == 0 && mycoords[ROW] == 0) { //I'm Proc. (0,0) 00122 for(int ip = 0; ip < dims[ROW]; ++ip){ 00123 int begidx = nbase[ip]*nx*ny; 00124 if (ip !=0){ // Proc (0,0) sends images and angle data to other processors 00125 srcoords[COL] = 0; 00126 srcoords[ROW] = ip; 00127 MPI_Cart_rank(comm_2d, srcoords, &otherpid); 00128 MPI_Send(&images[begidx],psize[ip]*nx*ny, MPI_FLOAT, otherpid, otherpid, comm_2d); 00129 MPI_Send(&angles[nbase[ip]*3],psize[ip]*3, MPI_FLOAT, otherpid, otherpid, comm_2d); 00130 } 00131 else{ // ip = 0: Proc (0,0) needs to copy images and angles into its imagesloc and anglesloc 00132 for (int i = 0; i < psize[ip]*nx*ny; i++){ 00133 imagesloc[i] = images[begidx+i]; 00134 } 00135 for (int i = 0; i < psize[ip]*3; i++){ 00136 anglesloc[i] = angles[nbase[ip]*3 + i]; 00137 } 00138 //printf("Finished copying to Proc (0,0) local"); 00139 } 00140 } //End for loop 00141 } //End if 00142 00143 if (mycoords[COL] == 0 && mycoords[ROW] != 0) { //I'm in g_c_0 and I'm not Processor (0,0) so I should receive data. 00144 MPI_Recv(imagesloc, psize[mycoords[ROW]]*nx*ny, MPI_FLOAT, 0, mypid, comm_2d, &mpistatus); 00145 MPI_Recv(anglesloc, psize[mycoords[ROW]]*3, MPI_FLOAT, 0, mypid, comm_2d, &mpistatus); 00146 } 00147 // Now have all the processors in group g_c_0 broadcast the images and angles along the row communicator 00148 srcoords[ROW] = 0; 00149 MPI_Cart_rank(comm_row, srcoords, &otherpid); 00150 MPI_Bcast(imagesloc, nangsloc*nx*ny, MPI_FLOAT, otherpid , comm_row); 00151 MPI_Bcast(anglesloc, nangsloc*3, MPI_FLOAT, otherpid , comm_row); 00152 00153 // Now distribute the volume (in spherical format) among columns of processors and use nnz to determine the splitting. Note: ptrs and coord are on all processors 00154 int radius; 00155 int volsize[3], origin[3]; 00156 volsize[0] = nx; 00157 volsize[1] = nx; 00158 volsize[2] = nx; 00159 origin[0] = nx/2+1; 00160 origin[1] = nx/2+1; 00161 origin[2] = nx/2+1; 00162 radius = nx/2-1; 00163 00164 ierr = getnnz( volsize, radius, origin, &nrays, &nnz); 00165 00166 int * ptrs = new int[nrays+1]; 00167 int * cord = new int[3*nrays]; 00168 ierr = getcb2sph(volsize, radius, origin, nnz, ptrs, cord); 00169 00170 int *nnzpart = new int[dims[COL]]; 00171 int *nnzbase = new int[dims[COL]+1]; 00172 nnzloc = setpart_gr1(comm_2d, nnz, nnzpart, nnzbase); 00173 00174 int *ptrstart = new int[dims[COL]+1]; 00175 nraysloc = sphpart(comm_2d, nrays, ptrs, nnzbase, ptrstart); 00176 00177 myptrstart = ptrstart[mycoords[COL]]; 00178 int nnzall[dims[COL]]; 00179 for (int i = 0; i<dims[COL]; i++) 00180 nnzall[i] = ptrs[ptrstart[i+1]] - ptrs[ptrstart[i]]; 00181 00182 nnzloc = nnzall[mycoords[COL]]; 00183 00184 // Print some stuff. 00185 printf("My coords are (%d,%d) and nangsloc = %d, nraysloc = %d, myptrstart = %d, nnzloc = %d\n", mycoords[ROW], mycoords[COL], nangsloc, nraysloc, myptrstart, nnzloc); 00186 00187 float *bvol_loc = new float[nnzloc]; 00188 float *vol_sphloc = new float[nnzloc]; 00189 for (int i=0; i< nnzloc; i++) 00190 bvol_loc[i] = 0.0; 00191 00192 // STEP 2: Have everyone perform the backprojection operation for their assigned images and portions of the volume. Then perform an Allreduce along the columns. 00193 00194 float phi, theta, psi; 00195 float dm[8]; 00196 00197 for (int i=0; i<nangsloc; i++){ 00198 phi = anglesloc[3*i+0]; 00199 theta = anglesloc[3*i+1]; 00200 psi = anglesloc[3*i+2]; 00201 dm[6] = 0; 00202 dm[7] = 0; 00203 00204 make_proj_mat(phi, theta, psi, dm); 00205 00206 ierr = bckpj3_Cart(volsize, nraysloc, nnzloc, dm, origin, radius, ptrs, cord, myptrstart, &imagesloc[nx*ny*i], bvol_loc); 00207 } 00208 00209 // Now an all reduce along the columns 00210 MPI_Allreduce (bvol_loc, vol_sphloc, nnzloc, MPI_FLOAT, MPI_SUM, comm_col); 00211 00212 // For testing purposes, we bring all the portions of the volume back together onto Proc (0,0). Note: we only need to deal with the first row of processors. 00213 if (mycoords[COL] != 0 && mycoords[ROW] == 0) { 00214 //Send data to Processor (0,0) 00215 srcoords[COL] = srcoords[ROW] = 0; 00216 MPI_Cart_rank(comm_2d, srcoords, &otherpid); 00217 MPI_Send(vol_sphloc, nnzloc, MPI_FLOAT, otherpid, otherpid, comm_2d); 00218 } 00219 float *onevol_sph = new float[nnz]; 00220 if (mycoords[COL] == 0 && mycoords[ROW] ==0){ 00221 //Copy data and recieve data 00222 float *vol_sph = new float[nnz]; 00223 for (int i=0; i<nnzloc; i++) 00224 vol_sph[i] = vol_sphloc[i]; 00225 00226 for (int i=1; i<dims[COL]; i++){ 00227 srcoords[ROW] = 0; 00228 srcoords[COL] = i; 00229 MPI_Cart_rank(comm_2d, srcoords, &otherpid); 00230 00231 MPI_Recv(&vol_sph[ptrs[ptrstart[i]]-1], nnzall[i], MPI_FLOAT, otherpid, mypid, comm_2d, &mpistatus); 00232 } 00233 //printf("Finished combining all volume parts\n"); 00234 00235 //Now compute the back projection serially on one processor (0,0) 00236 for (int i=0; i< nnz; i++) 00237 onevol_sph[i] = 0.0; 00238 00239 for (int i=0; i<nangs; i++){ 00240 phi = angles[3*i+0]; 00241 theta = angles[3*i+1]; 00242 psi = angles[3*i+2]; 00243 dm[6] = 0; 00244 dm[7] = 0; 00245 00246 make_proj_mat(phi, theta, psi, dm); 00247 00248 ierr = bckpj3(volsize, nrays, nnz, dm, origin, radius, ptrs, cord, &images[nx*ny*i], onevol_sph); 00249 } 00250 00251 float err=0; 00252 for (int i=0; i< nnz; i++){ 00253 err = err+(onevol_sph[i]-vol_sph[i])*(onevol_sph[i]-vol_sph[i]); 00254 } 00255 err = sqrt(err); 00256 printf("Cumulative error for backprojection is %f\n", err); 00257 delete [] vol_sph; 00258 } 00259 00260 // STEP 3: Now perform a forward projection operation for the assigned images and portions of the volume. Then perform an all_reduce along the rows. 00261 float * newimagesloc = new float[nangsloc*nx*ny]; 00262 for (int i=0; i<nangsloc*nx*ny; i++) 00263 newimagesloc[i] = 0.0; 00264 00265 for (int i=0; i<nangsloc; i++){ 00266 phi = anglesloc[3*i+0]; 00267 theta = anglesloc[3*i+1]; 00268 psi = anglesloc[3*i+2]; 00269 dm[6] = 0; 00270 dm[7] = 0; 00271 00272 make_proj_mat(phi, theta, psi, dm); 00273 00274 ierr = fwdpj3_Cart(volsize, nraysloc, nnzloc, dm, origin, radius, ptrs, cord, myptrstart, vol_sphloc, &newimagesloc[nx*ny*i]); 00275 } 00276 00277 // Now an all reduce along the rows 00278 MPI_Allreduce (newimagesloc, reprojloc, nangsloc*nx*ny, MPI_FLOAT, MPI_SUM, comm_row); 00279 00280 delete [] newimagesloc; 00281 // For testing purposes, we bring all the 2D images together onto Proc (0,0). Note: we only need to deal with the first column of processors. 00282 if (mycoords[ROW] != 0 && mycoords[COL] == 0) { 00283 //Send data to Processor (0,0) 00284 srcoords[COL] = srcoords[ROW] = 0; 00285 MPI_Cart_rank(comm_2d, srcoords, &otherpid); 00286 MPI_Send(reprojloc, nangsloc*nx*ny, MPI_FLOAT, otherpid, otherpid, comm_2d); 00287 } 00288 if (mycoords[COL] == 0 && mycoords[ROW] ==0){ 00289 //Copy data and recieve data 00290 float *reproj = new float[nangs*nx*ny]; 00291 for (int i=0; i<nangsloc*nx*ny; i++) 00292 reproj[i] = reprojloc[i]; 00293 00294 for (int i=1; i<dims[ROW]; i++){ 00295 srcoords[COL] = 0; 00296 srcoords[ROW] = i; 00297 MPI_Cart_rank(comm_2d, srcoords, &otherpid); 00298 00299 MPI_Recv(&reproj[nbase[i]*nx*ny], psize[i]*nx*ny, MPI_FLOAT, otherpid, mypid, comm_2d, &mpistatus); 00300 } 00301 delete [] reprojloc; 00302 // Now compute the forward projection serially on one processor (0,0) 00303 00304 float *allimages = new float[nangs*nx*ny]; 00305 for (int i=0; i< nangs*nx*ny; i++) 00306 allimages[i] = 0.0; 00307 00308 for (int i=0; i<nangs; i++){ 00309 phi = angles[3*i+0]; 00310 theta = angles[3*i+1]; 00311 psi = angles[3*i+2]; 00312 dm[6] = 0; 00313 dm[7] = 0; 00314 00315 make_proj_mat(phi, theta, psi, dm); 00316 00317 ierr = fwdpj3(volsize, nrays, nnz, dm, origin, radius, ptrs, cord, onevol_sph, &allimages[nx*ny*i]); 00318 00319 } 00320 00321 // Now compute the overall error. 00322 int idx; 00323 float err=0, max =0; 00324 for (int i=0; i< nangs*nx*ny; i++){ 00325 //if (allimages[i]!=reproj[i] && i < 256) 00326 // printf("i= %d\n",i); 00327 00328 err = err+(allimages[i]-reproj[i])*(allimages[i]-reproj[i]); 00329 00330 if (fabs(allimages[i]-reproj[i]) > max){ 00331 max = fabs(allimages[i]-reproj[i]); 00332 idx = i; 00333 } 00334 } 00335 err = sqrt(err); 00336 printf("Cumulative error for forward projection is %f with max error of %f occuring at %d\n", err, max, idx); 00337 printf("Max error: compare %f and %f\n",allimages[idx], reproj[idx]); 00338 00339 delete [] reproj; 00340 delete [] allimages; 00341 delete [] angles; 00342 delete [] images; 00343 } 00344 delete [] onevol_sph; 00345 delete [] vol_sphloc; 00346 delete [] bvol_loc; 00347 delete [] ptrs; 00348 delete [] cord; 00349 delete [] nnzpart; 00350 delete [] nnzbase; 00351 delete [] ptrstart; 00352 delete [] anglesloc; 00353 delete [] imagesloc; 00354 delete [] nbase; 00355 delete [] psize; 00356 00357 MPI_Comm_free(&comm_2d); 00358 MPI_Comm_free(&comm_row); 00359 MPI_Comm_free(&comm_col); 00360 00361 MPI_Finalize(); 00362 }
int make_proj_mat | ( | float | phi, | |
float | theta, | |||
float | psi, | |||
float * | dm | |||
) |
Definition at line 354 of file project3d.cpp.
00355 { 00356 // float cphi=cos(phi); 00357 // float sphi=sin(phi); 00358 // float cthe=cos(theta); 00359 // float sthe=sin(theta); 00360 // float cpsi=cos(psi); 00361 // float spsi=sin(psi); 00362 00363 double cphi, sphi, cthe, sthe, cpsi, spsi; 00364 double dphi = phi; 00365 double dthe = theta; 00366 double dpsi = psi; 00367 sincos(dphi, &sphi, &cphi); 00368 sincos(dthe, &sthe, &cthe); 00369 sincos(dpsi, &spsi, &cpsi); 00370 00371 dm[0]=cphi*cthe*cpsi-sphi*spsi; 00372 dm[1]=sphi*cthe*cpsi+cphi*spsi; 00373 dm[2]=-sthe*cpsi; 00374 dm[3]=-cphi*cthe*spsi-sphi*cpsi; 00375 dm[4]=-sphi*cthe*spsi+cphi*cpsi; 00376 dm[5]=sthe*spsi; 00377 00378 return 0; 00379 }
int setpart_gc1 | ( | MPI_Comm | comm_2d, | |
int | nangs, | |||
int * | psize, | |||
int * | nbase | |||
) |
Definition at line 323 of file utilcomm_Cart.cpp.
References COL, nangsloc, and ROW.
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 }
int setpart_gr1 | ( | MPI_Comm | comm_2d, | |
int | nnz, | |||
int * | nnzpart, | |||
int * | nnzbase | |||
) |
Definition at line 364 of file utilcomm_Cart.cpp.
References COL, nnzloc, and ROW.
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 }
int sphpart | ( | MPI_Comm | comm_2d, | |
int | nrays, | |||
int * | ptrs, | |||
int * | nnzbase, | |||
int * | ptrstart | |||
) |
Definition at line 85 of file project3d_Cart.cpp.
References COL, nraysloc, and ROW.
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 }
float* angles |
Definition at line 12 of file runcartrec.cpp.
float * anglesloc |
float * images |
Definition at line 10 of file runcartrec.cpp.
float* imagesloc |
int myptrstart |
Definition at line 8 of file runcartrec.cpp.
Referenced by main(), recons3d_CGLS_mpi_Cart(), recons3d_HyBR_mpi_Cart(), and recons3d_sirt_mpi_Cart().
int nangs |
int nangsloc |
float* newimagesloc |
int nnz |
Definition at line 9 of file runcartrec.cpp.
Referenced by ali3d_d(), EMAN::ChaoProjector::backproject3d(), EMAN::ChaoProjector::cb2sph(), cb2sph(), wustl_mm::SkeletonMaker::Volume::countInt(), wustl_mm::SkeletonMaker::Volume::countIntEuler(), getcb2sph(), EMAN::ChaoProjector::getnnz(), main(), EMAN::ChaoProjector::project3d(), recons3d_CGLS_mpi_Cart(), recons3d_HyBR_mpi_Cart(), recons3d_sirt_mpi(), recons3d_sirt_mpi_Cart(), EMAN::ChaoProjector::sph2cb(), sph2cb(), and unified().
int nnzloc |
Definition at line 9 of file runcartrec.cpp.
Referenced by main(), recons3d_CGLS_mpi_Cart(), recons3d_HyBR_mpi_Cart(), recons3d_sirt_mpi_Cart(), and setpart_gr1().
int nrays |
Definition at line 9 of file runcartrec.cpp.
Referenced by ali3d_d(), EMAN::ChaoProjector::backproject3d(), EMAN::ChaoProjector::cb2sph(), cb2sph(), getcb2sph(), EMAN::ChaoProjector::getnnz(), main(), EMAN::ChaoProjector::project3d(), recons3d_CGLS_mpi_Cart(), recons3d_HyBR_mpi_Cart(), recons3d_sirt_mpi(), recons3d_sirt_mpi_Cart(), and unified().
int nraysloc |
Definition at line 9 of file runcartrec.cpp.
Referenced by main(), recons3d_CGLS_mpi_Cart(), recons3d_HyBR_mpi_Cart(), recons3d_sirt_mpi_Cart(), and sphpart().
int nx |
Definition at line 8 of file runcartrec.cpp.
Referenced by add_complex_at_fast(), EMAN::newfile_store::add_image(), EMAN::FourierWeightAverager::add_image(), EMAN::ImageAverager::add_image(), EMAN::TomoAverager::add_image(), EMAN::Util::add_img(), EMAN::Util::add_img2(), EMAN::Util::add_img_abs(), EMAN::Util::addn_img(), ali3d_d(), EMAN::FRM2DAligner::align(), EMAN::RTFSlowExhaustiveAligner::align(), EMAN::TranslationalAligner::align(), EMAN::PCAlarge::analyze(), apmd(), apmq(), EMAN::BoxingTools::auto_correlation_pick(), EMAN::PawelProjector::backproject3d(), EMAN::ChaoProjector::bckpj3(), bckpj3(), bckpj3_Cart(), EMAN::ChaoProjector::cb2sph(), cb2sph(), CleanStack(), CleanStack_Cart(), EMAN::Util::cluster_equalsize(), EMAN::Util::cluster_pairwise(), EMAN::Util::cml_line_in3d(), EMAN::Util::cml_prepare_line(), EMAN::XYZCmp::cmp(), EMAN::OptVarianceCmp::cmp(), EMAN::QuadMinDotCmp::cmp(), EMAN::TomoFscCmp::cmp(), EMAN::DotCmp::cmp(), EMAN::SqEuclideanCmp::cmp(), cmplx(), wustl_mm::SkeletonMaker::Volume::components26(), wustl_mm::SkeletonMaker::Volume::components6(), EMAN::Util::compress_image_mask(), EMAN::EMAN2Ctf::compute_2d_complex(), EMAN::EMAN1Ctf::compute_2d_complex(), wustl_mm::SkeletonMaker::Volume::countInt(), wustl_mm::SkeletonMaker::Volume::countIntEuler(), wustl_mm::SkeletonMaker::Volume::curveSkeleton(), wustl_mm::SkeletonMaker::Volume::curveSkeleton2D(), EMAN::Util::cyclicshift(), EMAN::Util::decimate(), EMAN::Util::div_filter(), EMAN::Util::div_img(), EMAN::Util::divn_filter(), EMAN::Util::divn_img(), EMAN::PCA::dopca_lan(), EMAN::PCA::dopca_ooc(), EMAN::TestUtil::dump_emdata(), EMAN::Processor::EMFourierFilterFunc(), wustl_mm::SkeletonMaker::Volume::erodeHelix(), wustl_mm::SkeletonMaker::Volume::erodeSheet(), fcalc(), fgcalc(), find_group(), EMAN::FourierWeightAverager::finish(), EMAN::TomoAverager::finish(), EMAN::fourierproduct(), EMAN::ChaoProjector::fwdpj3(), fwdpj3(), fwdpj3_Cart(), EMAN::Util::get_biggest_cluster(), get_complex_index_fast(), get_size(), EMAN::Util::get_slice(), get_value_at(), get_xsize(), EMAN::XYData::get_yatx(), getcb2sph(), EMAN::ChaoProjector::getnnz(), getnnz(), wustl_mm::SkeletonMaker::Volume::getNumNeighbor6(), wustl_mm::SkeletonMaker::Volume::getNumPotComplex(), wustl_mm::SkeletonMaker::Volume::hasCompleteHelix(), wustl_mm::SkeletonMaker::Volume::hasCompleteSheet(), EMAN::Util::histogram(), ilaenv_(), EMAN::Util::im_diff(), EMAN::Util::image_mutation(), EMAN::Util::infomask(), EMAN::ImagicIO2::init_test(), EMAN::PCAlarge::insert_image(), EMAN::PCAsmall::insert_image(), EMAN::GaussFFTProjector::interp_ft_3d(), EMAN::OmapIO::is_valid(), EMAN::MrcIO::is_valid(), EMAN::ImagicIO2::is_valid(), EMAN::ImagicIO::is_valid(), EMAN::EmIO::is_valid(), wustl_mm::SkeletonMaker::Volume::isFeatureFace(), wustl_mm::SkeletonMaker::Volume::isHelixEnd(), wustl_mm::SkeletonMaker::Volume::isSheetEnd(), EMAN::Util::mad_scalar(), EMAN::Util::madn_scalar(), main(), EMAN::TestUtil::make_image_file(), EMAN::TestUtil::make_image_file2(), wustl_mm::SkeletonMaker::Volume::markCellFace(), EMAN::Util::min_dist_four(), EMAN::Util::move_points(), EMAN::Util::mul_img(), EMAN::Util::mul_scalar(), EMAN::Util::muln_img(), EMAN::Util::mult_scalar(), operator()(), EMAN::Util::pack_complex_to_real(), EMAN::Util::pad(), EMAN::padfft_slice(), EMAN::periodogram(), EMAN::Util::Polar2Dmi(), printImage(), EMAN::Util::printMatI3D(), EMAN::KmeansSegmentProcessor::process(), EMAN::DistanceSegmentProcessor::process(), EMAN::CutoffBlockProcessor::process_inplace(), EMAN::DiffBlockProcessor::process_inplace(), EMAN::BoxStatProcessor::process_inplace(), EMAN::PaintProcessor::process_inplace(), EMAN::LowpassRandomPhaseProcessor::process_inplace(), EMAN::LinearPyramidProcessor::process_inplace(), EMAN::FourierGriddingProjector::project3d(), EMAN::MaxValProjector::project3d(), EMAN::StandardProjector::project3d(), EMAN::PawelProjector::project3d(), ReadStackandDist(), ReadStackandDist_Cart(), ReadVandBcast(), recons3d_CGLS_mpi_Cart(), recons3d_HyBR_mpi_Cart(), recons3d_sirt_mpi(), recons3d_sirt_mpi_Cart(), EMAN::Util::reconstitute_image_mask(), set_data(), EMAN::PointArray::set_from_density_map(), EMAN::Util::set_line(), EMAN::varimax::set_params(), EMAN::PCAlarge::set_params(), EMAN::PCAsmall::set_params(), set_value_at(), set_value_at_fast(), sgebrd_(), sgelqf_(), sgeqrf_(), wustl_mm::SkeletonMaker::Volume::skeleton(), sorglq_(), sorgql_(), sorgqr_(), EMAN::ChaoProjector::sph2cb(), sph2cb(), ssytrd_(), EMAN::Util::sub_img(), EMAN::Util::subn_img(), wustl_mm::SkeletonMaker::Volume::surfaceSkeletonPres(), EMAN::TestUtil::to_emobject(), unified(), EMAN::TestUtil::verify_image_file(), EMAN::TestUtil::verify_image_file2(), EMAN::EMUtil::vertical_acf(), EMAN::Util::window(), and EMAN::SpiderIO::write_single_header().
int ny |
Definition at line 8 of file runcartrec.cpp.
Referenced by add_complex_at_fast(), EMAN::newfile_store::add_image(), EMAN::FourierWeightAverager::add_image(), EMAN::ImageAverager::add_image(), EMAN::TomoAverager::add_image(), EMAN::Util::add_img(), EMAN::Util::add_img2(), EMAN::Util::add_img_abs(), EMAN::Util::addn_img(), ali3d_d(), EMAN::RTFExhaustiveAligner::align(), EMAN::RotatePrecenterAligner::align(), EMAN::TranslationalAligner::align(), apmd(), apmq(), EMAN::BoxingTools::auto_correlation_pick(), EMAN::PawelProjector::backproject3d(), bckpj3(), bckpj3_Cart(), EMAN::ChaoProjector::cb2sph(), cb2sph(), CleanStack(), CleanStack_Cart(), EMAN::Util::cml_line_in3d(), EMAN::XYZCmp::cmp(), EMAN::FRCCmp::cmp(), EMAN::PhaseCmp::cmp(), EMAN::QuadMinDotCmp::cmp(), EMAN::TomoFscCmp::cmp(), EMAN::DotCmp::cmp(), EMAN::SqEuclideanCmp::cmp(), cmplx(), wustl_mm::SkeletonMaker::Volume::components26(), wustl_mm::SkeletonMaker::Volume::components6(), EMAN::Util::compress_image_mask(), EMAN::EMAN2Ctf::compute_2d_complex(), EMAN::EMAN1Ctf::compute_2d_complex(), wustl_mm::SkeletonMaker::Volume::countInt(), wustl_mm::SkeletonMaker::Volume::countIntEuler(), wustl_mm::SkeletonMaker::Volume::curveSkeleton(), wustl_mm::SkeletonMaker::Volume::curveSkeleton2D(), EMAN::Util::cyclicshift(), EMAN::Util::decimate(), EMAN::Util::div_filter(), EMAN::Util::div_img(), EMAN::Util::divn_filter(), EMAN::Util::divn_img(), EMAN::TestUtil::dump_emdata(), EMAN::Processor::EMFourierFilterFunc(), wustl_mm::SkeletonMaker::Volume::erodeHelix(), wustl_mm::SkeletonMaker::Volume::erodeSheet(), fcalc(), fgcalc(), find_group(), EMAN::FourierWeightAverager::finish(), EMAN::TomoAverager::finish(), EMAN::fourierproduct(), fwdpj3(), fwdpj3_Cart(), EMAN::Util::get_biggest_cluster(), get_complex_index_fast(), get_ndim(), get_size(), EMAN::Util::get_slice(), get_ysize(), getcb2sph(), EMAN::ChaoProjector::getnnz(), getnnz(), wustl_mm::SkeletonMaker::Volume::getNumNeighbor6(), wustl_mm::SkeletonMaker::Volume::getNumPotComplex(), wustl_mm::SkeletonMaker::Volume::hasCompleteHelix(), wustl_mm::SkeletonMaker::Volume::hasCompleteSheet(), EMAN::Util::histogram(), EMAN::Util::im_diff(), EMAN::Util::infomask(), EMAN::ImagicIO2::init_test(), EMAN::GaussFFTProjector::interp_ft_3d(), EMAN::OmapIO::is_valid(), EMAN::MrcIO::is_valid(), EMAN::ImagicIO2::is_valid(), EMAN::ImagicIO::is_valid(), EMAN::EmIO::is_valid(), wustl_mm::SkeletonMaker::Volume::isFeatureFace(), wustl_mm::SkeletonMaker::Volume::isHelixEnd(), wustl_mm::SkeletonMaker::Volume::isSheetEnd(), EMAN::Util::mad_scalar(), EMAN::Util::madn_scalar(), main(), EMAN::TestUtil::make_image_file(), EMAN::TestUtil::make_image_file2(), wustl_mm::SkeletonMaker::Volume::markCellFace(), EMAN::Util::min_dist_four(), EMAN::Util::move_points(), EMAN::Util::mul_img(), EMAN::Util::mul_scalar(), EMAN::Util::muln_img(), EMAN::Util::mult_scalar(), operator()(), EMAN::Util::pack_complex_to_real(), EMAN::Util::pad(), EMAN::padfft_slice(), EMAN::periodogram(), EMAN::Util::Polar2Dmi(), printImage(), EMAN::Util::printMatI3D(), EMAN::KmeansSegmentProcessor::process(), EMAN::DistanceSegmentProcessor::process(), EMAN::CutoffBlockProcessor::process_inplace(), EMAN::DiffBlockProcessor::process_inplace(), EMAN::BoxStatProcessor::process_inplace(), EMAN::PaintProcessor::process_inplace(), EMAN::LowpassRandomPhaseProcessor::process_inplace(), EMAN::LinearPyramidProcessor::process_inplace(), EMAN::FourierGriddingProjector::project3d(), EMAN::MaxValProjector::project3d(), EMAN::StandardProjector::project3d(), EMAN::PawelProjector::project3d(), ReadStackandDist(), ReadStackandDist_Cart(), ReadVandBcast(), EMAN::Util::reconstitute_image_mask(), set_complex_size(), set_data(), EMAN::PointArray::set_from_density_map(), EMAN::varimax::set_params(), EMAN::PCAlarge::set_params(), EMAN::PCAsmall::set_params(), set_value_at(), wustl_mm::SkeletonMaker::Volume::skeleton(), EMAN::ChaoProjector::sph2cb(), sph2cb(), EMAN::Util::sub_img(), EMAN::Util::subn_img(), wustl_mm::SkeletonMaker::Volume::surfaceSkeletonPres(), EMAN::TestUtil::to_emobject(), unified(), EMAN::Util::vareas(), EMAN::TestUtil::verify_image_file(), EMAN::TestUtil::verify_image_file2(), EMAN::EMUtil::vertical_acf(), EMAN::Util::window(), and EMAN::SpiderIO::write_single_header().
float* onevol_sph |
float * reproj |
float * reprojloc |
float* vol_sph |