#include <stdio.h>
#include <stdlib.h>
#include <mpi.h>
#include <math.h>
#include <string.h>
Include dependency graph for utilcomm2d.cpp:
Go to the source code of this file.
Defines | |
#define | cord(i, j) cord[((j)-1)*3 + (i) -1] |
#define | ptrs(i) ptrs[(i)-1] |
#define | dm(i) dm[(i)-1] |
#define | x(i) x[(i)-1] |
#define | y(i, j) y[((j)-1)*nx + (i) - 1] |
#define | y(i) y[(i)-1] |
#define | x(i, j) x[((j)-1)*nx + (i) - 1] |
#define | x(i) x[(i)-1] |
#define | y(i, j) y[((j)-1)*nx + (i) - 1] |
#define | y(i) y[(i)-1] |
#define | x(i, j) x[((j)-1)*nx + (i) - 1] |
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) |
#define cord | ( | i, | |||
j | ) | cord[((j)-1)*3 + (i) -1] |
Definition at line 140 of file utilcomm2d.cpp.
#define dm | ( | i | ) | dm[(i)-1] |
Definition at line 142 of file utilcomm2d.cpp.
#define ptrs | ( | i | ) | ptrs[(i)-1] |
Definition at line 141 of file utilcomm2d.cpp.
#define x | ( | i, | |||
j | ) | x[((j)-1)*nx + (i) - 1] |
Definition at line 615 of file utilcomm2d.cpp.
#define x | ( | i | ) | x[(i)-1] |
Definition at line 615 of file utilcomm2d.cpp.
#define x | ( | i, | |||
j | ) | x[((j)-1)*nx + (i) - 1] |
Definition at line 615 of file utilcomm2d.cpp.
#define x | ( | i | ) | x[(i)-1] |
Definition at line 615 of file utilcomm2d.cpp.
#define y | ( | i | ) | y[(i)-1] |
Definition at line 614 of file utilcomm2d.cpp.
#define y | ( | i, | |||
j | ) | y[((j)-1)*nx + (i) - 1] |
Definition at line 614 of file utilcomm2d.cpp.
#define y | ( | i | ) | y[(i)-1] |
Definition at line 614 of file utilcomm2d.cpp.
#define y | ( | i, | |||
j | ) | y[((j)-1)*nx + (i) - 1] |
Definition at line 614 of file utilcomm2d.cpp.
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 305 of file utilcomm2d.cpp.
00308 { 00309 int ia = 0; 00310 if (a >= 0.0) { 00311 ia = (int)floor(a); 00312 } 00313 else { 00314 ia = (int)ceil(a); 00315 } 00316 return ia; 00317 }
int make_proj_mat | ( | float | phi, | |
float | theta, | |||
float | psi, | |||
float * | dm | |||
) |
Definition at line 272 of file utilcomm2d.cpp.
00274 : 00275 phi, theta, psi - Euler angles 00276 Output: 00277 dm - projection matrix 00278 */ 00279 { 00280 // float cphi=cos(phi); 00281 // float sphi=sin(phi); 00282 // float cthe=cos(theta); 00283 // float sthe=sin(theta); 00284 // float cpsi=cos(psi); 00285 // float spsi=sin(psi); 00286 00287 double cphi, sphi, cthe, sthe, cpsi, spsi; 00288 double dphi = phi; 00289 double dthe = theta; 00290 double dpsi = psi; 00291 sincos(dphi, &sphi, &cphi); 00292 sincos(dthe, &sthe, &cthe); 00293 sincos(dpsi, &spsi, &cpsi); 00294 00295 dm[0]=cphi*cthe*cpsi-sphi*spsi; 00296 dm[1]=sphi*cthe*cpsi+cphi*spsi; 00297 dm[2]=-sthe*cpsi; 00298 dm[3]=-cphi*cthe*spsi-sphi*cpsi; 00299 dm[4]=-sphi*cthe*spsi+cphi*cpsi; 00300 dm[5]=sthe*spsi; 00301 00302 return 0; 00303 }
int setpart_gc1 | ( | MPI_Comm | comm_2d, | |
int | nangs, | |||
int * | psize, | |||
int * | nbase | |||
) |
Definition at line 7 of file utilcomm2d.cpp.
References COL, nangsloc, and ROW.
00010 : 00011 comm_2d - Cartesian communicator 00012 nangs - number of 2D images 00013 Output: 00014 psize - vector containing the partition distribution 00015 nbase - vector containing the base of the partitions 00016 nangsloc - number of local angles 00017 */ 00018 { 00019 int ROW = 0, COL = 1; 00020 int dims[2], periods[2], mycoords[2]; 00021 int nangsloc, nrem; 00022 00023 // Get information associated with comm_2d 00024 MPI_Cart_get(comm_2d, 2, dims, periods, mycoords); 00025 00026 nangsloc = nangs/dims[ROW]; 00027 nrem = nangs - dims[ROW]*nangsloc; 00028 if (mycoords[ROW] < nrem) nangsloc++; 00029 00030 for (int i = 0; i < dims[ROW]; i++) { 00031 psize[i] = nangs/dims[ROW]; 00032 if (i < nrem) psize[i]++; 00033 } 00034 00035 nbase[0] = 0; 00036 for (int i = 1; i < dims[ROW]; i++) { 00037 nbase[i] = nbase[i-1] + psize[i-1]; 00038 } 00039 00040 return nangsloc; 00041 }
int setpart_gr1 | ( | MPI_Comm | comm_2d, | |
int | nnz, | |||
int * | nnzpart, | |||
int * | nnzbase | |||
) |
Definition at line 104 of file utilcomm2d.cpp.
References COL, nnzloc, and ROW.
00107 : 00108 comm_2d - Cartesian communicator 00109 nnz - total number of non-zeros in the volume 00110 Output: 00111 nnzpart - vector containing the partition distribution 00112 nnzbase - vector containing the base of the partitions 00113 nnzloc - initial local nnz 00114 */ 00115 { 00116 int ROW = 0, COL = 1; 00117 int dims[2], periods[2], mycoords[2]; 00118 int nnzloc, nrem; 00119 00120 // Get information associated with comm_2d 00121 MPI_Cart_get(comm_2d, 2, dims, periods, mycoords); 00122 00123 nnzloc = nnz/dims[COL]; 00124 nrem = nnz - dims[COL]*nnzloc; 00125 if (mycoords[COL] < nrem) nnzloc++; 00126 00127 for (int i = 0; i < dims[COL]; i++) { 00128 nnzpart[i] = nnz/dims[COL]; 00129 if (i < nrem) nnzpart[i]++; 00130 } 00131 00132 nnzbase[0] = 0; 00133 for (int i = 1; i < dims[COL]; i++) { 00134 nnzbase[i] = nnzbase[i-1] + nnzpart[i-1]; 00135 } 00136 nnzbase[dims[COL]] = nnz; 00137 return nnzloc; 00138 }
int sphpart | ( | MPI_Comm | comm_2d, | |
int | nrays, | |||
int * | ptrs, | |||
int * | nnzbase, | |||
int * | ptrstart | |||
) |
Definition at line 214 of file utilcomm2d.cpp.
References COL, nraysloc, and ROW.
00217 : 00218 comm_2d - Cartesian communicator 00219 nrays - total number of rays 00220 ptrs - vector containing pointers 00221 nnzbase - ideal volume partition of nnz 00222 Output: 00223 ptrstart - vector containing all the starting pointers for each column processor group 00224 nraysloc - actual number of local rays 00225 */ 00226 { 00227 int ROW = 0, COL = 1; 00228 int dims[2], periods[2], mycoords[2]; 00229 int nraysloc = 0; 00230 00231 // Get information associated with comm_2d 00232 MPI_Cart_get(comm_2d, 2, dims, periods, mycoords); 00233 00234 int count = 1; 00235 if (mycoords[COL] == 0){ // First column starts out with the first ray 00236 nraysloc++; 00237 } 00238 ptrstart[0] = 0; 00239 00240 for(int i=1; i<nrays ; i++){ 00241 if (ptrs[i]-1 <= nnzbase[count] && ptrs[i+1]-1 >= nnzbase[count]){ 00242 //nnzbase is between or equal to ptrs 00243 00244 if (nnzbase[count]-(ptrs[i]-1)>= ptrs[i+1]-1-nnzbase[count]){ 00245 if(mycoords[COL] == count-1){ // ray belongs to count-1 00246 nraysloc++; 00247 } 00248 ptrstart[count] = i+1; 00249 count++; 00250 00251 } else { //nnzbase[count]-(ptrs[i]-1)>= ptrs[i+1]-1-nnzbase[count] 00252 if(mycoords[COL] == count){// ray belongs to count and it's a new starting ray 00253 nraysloc++; 00254 } 00255 ptrstart[count] = i; 00256 count++; 00257 } 00258 00259 } 00260 else { //ptrs[i]-1 > nnzbase[count] so ray belongs to count-1 00261 if(mycoords[COL] == count-1){ 00262 nraysloc++; 00263 } 00264 00265 } 00266 } // end for 00267 ptrstart[dims[COL]] = nrays; 00268 return nraysloc; 00269 00270 }