#include <projector.h>
Inheritance diagram for EMAN::ChaoProjector:
Public Member Functions | |
EMData * | project3d (EMData *vol) const |
Project an 3D image into a 2D image. | |
EMData * | backproject3d (EMData *imagestack) const |
Back-project a 2D image into a 3D image. | |
string | get_name () const |
Get the projector's name. | |
string | get_desc () const |
TypeDict | get_param_types () const |
Get processor parameter information in a dictionary. | |
Static Public Member Functions | |
static Projector * | NEW () |
Static Public Attributes | |
static const string | NAME = "chao" |
Private Member Functions | |
int | getnnz (Vec3i volsize, int ri, Vec3i origin, int *nray, int *nnz) const |
int | cb2sph (float *cube, Vec3i volsize, int ri, Vec3i origin, int nnz0, int *ptrs, int *cord, float *sphere) const |
int | sph2cb (float *sphere, Vec3i volsize, int nray, int ri, int nnz0, int *ptrs, int *cord, float *cube) const |
int | fwdpj3 (Vec3i volsize, int nray, int nnz, float *dm, Vec3i origin, int ri, int *ptrs, int *cord, float *x, float *y) const |
int | bckpj3 (Vec3i volsize, int nray, int nnz, float *dm, Vec3i origin, int ri, int *ptrs, int *cord, float *x, float *y) const |
int | ifix (float a) const |
void | setdm (vector< float > anglelist, string const angletype, float *dm) const |
(C. Yang)
Definition at line 338 of file projector.h.
Back-project a 2D image into a 3D image.
Implements EMAN::Projector.
Definition at line 1721 of file projector.cpp.
References anglelist, bckpj3(), cb2sph(), cord, cube, dm, EMDeleteArray(), EMAN::EMData::get_data(), EMAN::Util::get_min(), EMAN::Transform::get_rotation(), EMAN::EMData::get_xsize(), EMAN::EMData::get_ysize(), EMAN::EMData::get_zsize(), getnnz(), EMAN::Dict::has_key(), images, LOGERR, nnz, nrays, NullPointerException, EMAN::Projector::params, phi, ptrs, EMAN::EMData::set_complex(), EMAN::EMData::set_ri(), EMAN::EMData::set_size(), setdm(), EMAN::Dict::size(), sph2cb(), sphere, status, theta, EMAN::EMData::to_zero(), and EMAN::EMData::update().
01722 { 01723 int nrays, nnz, status, j; 01724 float *dm; 01725 int *ptrs, *cord; 01726 float *sphere, *images, *cube; 01727 01728 int nximg = imagestack->get_xsize(); 01729 int nyimg = imagestack->get_ysize(); 01730 int nslices = imagestack->get_zsize(); 01731 01732 int dim = Util::get_min(nximg,nyimg); 01733 Vec3i volsize(nximg,nyimg,dim); 01734 01735 Vec3i origin(0,0,0); 01736 // If a sensible origin isn't passed in, choose the middle of 01737 // the cube. 01738 if (params.has_key("origin_x")) {origin[0] = params["origin_x"];} 01739 else {origin[0] = nximg/2+1;} 01740 if (params.has_key("origin_y")) {origin[1] = params["origin_y"];} 01741 else {origin[1] = nyimg/2+1;} 01742 if (params.has_key("origin_z")) {origin[1] = params["origin_z"];} 01743 else {origin[2] = dim/2+1;} 01744 01745 int ri; 01746 if (params.has_key("radius")) {ri = params["radius"];} 01747 else {ri = dim/2 - 1;} 01748 01749 // retrieve the voxel values 01750 images = imagestack->get_data(); 01751 01752 // count the number of voxels within a sphere centered at icent, 01753 // with radius ri 01754 status = getnnz(volsize, ri, origin, &nrays, &nnz); 01755 // need to check status... 01756 01757 // convert from cube to sphere 01758 sphere = new float[nnz]; 01759 ptrs = new int[nrays+1]; 01760 cord = new int[3*nrays]; 01761 if (sphere == NULL || ptrs == NULL || cord == NULL) { 01762 fprintf(stderr,"ChaoProjector::backproject3d, failed to allocate!\n"); 01763 exit(1); 01764 } 01765 for (int i = 0; i<nnz; i++) sphere[i] = 0.0; 01766 for (int i = 0; i<nrays+1; i++) ptrs[i] = 0; 01767 for (int i = 0; i<3*nrays; i++) cord[i] = 0; 01768 01769 int nangles = 0; 01770 vector<float> anglelist; 01771 string angletype = "SPIDER"; 01772 // Do we have a list of angles? 01773 if (params.has_key("anglelist")) { 01774 anglelist = params["anglelist"]; 01775 nangles = anglelist.size() / 3; 01776 } else { 01777 Transform* t3d = params["transform"]; 01778 if ( t3d == NULL ) throw NullPointerException("The transform object (required for projection), was not specified"); 01779 // This part was modified by David Woolford - 01780 // Before this the code worked only for SPIDER and EMAN angles, 01781 // but the framework of the Transform3D allows for a generic implementation 01782 // as specified here. 01783 // This was broken by david. we need here a loop over all projections and put all angles on stack PAP 06/28/09 01784 Dict p = t3d->get_rotation("spider"); 01785 if(t3d) {delete t3d; t3d=0;} 01786 01787 float phi = p["phi"]; 01788 float theta = p["theta"]; 01789 float psi = p["psi"]; 01790 anglelist.push_back(phi); 01791 anglelist.push_back(theta); 01792 anglelist.push_back(psi); 01793 nangles = 1; 01794 } 01795 01796 // End David Woolford modifications 01797 01798 if (nslices != nangles) { 01799 LOGERR("the number of images does not match the number of angles"); 01800 return 0; 01801 } 01802 01803 dm = new float[nangles*9]; 01804 setdm(anglelist, angletype, dm); 01805 01806 // return volume 01807 EMData *ret = new EMData(); 01808 ret->set_size(nximg, nyimg, dim); 01809 ret->set_complex(false); 01810 ret->set_ri(true); 01811 ret->to_zero(); 01812 01813 cube = ret->get_data(); 01814 // cb2sph should be replaced by something that touches only ptrs and cord 01815 status = cb2sph(cube, volsize, ri, origin, nnz, ptrs, cord, sphere); 01816 // check status 01817 01818 for (j = 1; j <= nangles; j++) { 01819 status = bckpj3(volsize, nrays, nnz, &dm(1,j), origin, ri, 01820 ptrs , cord , &images(1,1,j), sphere); 01821 // check status? 01822 } 01823 01824 status = sph2cb(sphere, volsize, nrays, ri, nnz, ptrs, cord, cube); 01825 // check status? 01826 01827 // deallocate all temporary work space 01828 EMDeleteArray(dm); 01829 EMDeleteArray(ptrs); 01830 EMDeleteArray(cord); 01831 EMDeleteArray(sphere); 01832 01833 ret->update(); 01834 return ret; 01835 }
int ChaoProjector::bckpj3 | ( | Vec3i | volsize, | |
int | nray, | |||
int | nnz, | |||
float * | dm, | |||
Vec3i | origin, | |||
int | ri, | |||
int * | ptrs, | |||
int * | cord, | |||
float * | x, | |||
float * | y | |||
) | const [private] |
Definition at line 1474 of file projector.cpp.
References cord, dm, ifix(), nx, ptrs, status, x, and y.
Referenced by backproject3d().
01477 { 01478 int i, j, iqx,iqy, xc, yc, zc; 01479 float xb, yb, dx, dy, dx1m, dy1m, dxdy; 01480 int status = 0; 01481 01482 int xcent = origin[0]; 01483 int ycent = origin[1]; 01484 int zcent = origin[2]; 01485 01486 int nx = volsize[0]; 01487 01488 if ( nx > 2*ri) { 01489 for (i = 1; i <= nrays; i++) { 01490 zc = cord(1,i) - zcent; 01491 yc = cord(2,i) - ycent; 01492 xc = cord(3,i) - xcent; 01493 01494 xb = zc*dm(1)+yc*dm(2)+xc*dm(3) + xcent; 01495 yb = zc*dm(4)+yc*dm(5)+xc*dm(6) + ycent; 01496 01497 for (j = ptrs(i); j <ptrs(i+1); j++) { 01498 iqx = ifix((float)(xb)); 01499 iqy = ifix((float)(yb)); 01500 01501 dx = xb - (float)(iqx); 01502 dy = yb - (float)(iqy); 01503 dx1m = 1.0f - dx; 01504 dy1m = 1.0f - dy; 01505 dxdy = dx*dy; 01506 /* 01507 c y(j) = y(j) + dx1m*dy1m*x(iqx , iqy) 01508 c & + dx1m*dy *x(iqx , iqy+1) 01509 c & + dx *dy1m*x(iqx+1, iqy) 01510 c & + dx *dy *x(iqx+1, iqy+1) 01511 c 01512 c --- faster version of the above commented out 01513 c code (derived by summing the following table 01514 c of coefficients along the colunms) --- 01515 c 01516 c 1 dx dy dxdy 01517 c ------ -------- -------- ------- 01518 c x(i,j) -x(i,j) -x(i,j) x(i,j) 01519 c x(i,j+1) -x(i,j+1) 01520 c x(i+1,j) -x(i+1,j) 01521 c x(i+1,j+1) 01522 c 01523 */ 01524 y(j) += x(iqx,iqy) 01525 + dx*(-x(iqx,iqy)+x(iqx+1,iqy)) 01526 + dy*(-x(iqx,iqy)+x(iqx,iqy+1)) 01527 + dxdy*( x(iqx,iqy) - x(iqx,iqy+1) 01528 -x(iqx+1,iqy) + x(iqx+1,iqy+1) ); 01529 01530 xb += dm(1); 01531 yb += dm(4); 01532 } // end for j 01533 } // end for i 01534 } 01535 else { 01536 fprintf(stderr, "bckpj3: nx must be greater than 2*ri\n"); 01537 } 01538 01539 return status; 01540 }
int ChaoProjector::cb2sph | ( | float * | cube, | |
Vec3i | volsize, | |||
int | ri, | |||
Vec3i | origin, | |||
int | nnz0, | |||
int * | ptrs, | |||
int * | cord, | |||
float * | sphere | |||
) | const [private] |
Definition at line 1304 of file projector.cpp.
References cord, cube, nnz, nrays, nx, ny, ptrs, sphere, and status.
Referenced by backproject3d(), and project3d().
01306 { 01307 int xs, ys, zs, xx, yy, zz, rs, r2; 01308 int ix, iy, iz, jnz, nnz, nrays; 01309 int ftm = 0, status = 0; 01310 01311 int xcent = (int)origin[0]; 01312 int ycent = (int)origin[1]; 01313 int zcent = (int)origin[2]; 01314 01315 int nx = (int)volsize[0]; 01316 int ny = (int)volsize[1]; 01317 int nz = (int)volsize[2]; 01318 01319 r2 = ri*ri; 01320 nnz = 0; 01321 nrays = 0; 01322 ptrs(1) = 1; 01323 01324 for (ix = 1; ix <= nx; ix++) { 01325 xs = ix-xcent; 01326 xx = xs*xs; 01327 for ( iy = 1; iy <= ny; iy++ ) { 01328 ys = iy-ycent; 01329 yy = ys*ys; 01330 jnz = 0; 01331 01332 ftm = 1; 01333 // not the most efficient implementation 01334 for (iz = 1; iz <= nz; iz++) { 01335 zs = iz-zcent; 01336 zz = zs*zs; 01337 rs = xx + yy + zz; 01338 if (rs <= r2) { 01339 jnz++; 01340 nnz++; 01341 sphere(nnz) = cube(iz, iy, ix); 01342 01343 // record the coordinates of the first nonzero === 01344 if (ftm) { 01345 nrays++; 01346 cord(1,nrays) = iz; 01347 cord(2,nrays) = iy; 01348 cord(3,nrays) = ix; 01349 ftm = 0; 01350 } 01351 } 01352 } // end for (iz..) 01353 if (jnz > 0) { 01354 ptrs(nrays+1) = ptrs(nrays) + jnz; 01355 } // endif (jnz) 01356 } // end for iy 01357 } // end for ix 01358 if (nnz != nnz0) status = -1; 01359 return status; 01360 }
int ChaoProjector::fwdpj3 | ( | Vec3i | volsize, | |
int | nray, | |||
int | nnz, | |||
float * | dm, | |||
Vec3i | origin, | |||
int | ri, | |||
int * | ptrs, | |||
int * | cord, | |||
float * | x, | |||
float * | y | |||
) | const [private] |
Definition at line 1398 of file projector.cpp.
References cord, dm, ifix(), nx, ptrs, status, x, and y.
Referenced by project3d().
01401 { 01402 /* 01403 purpose: y <--- proj(x) 01404 input : volsize the size (nx,ny,nz) of the volume 01405 nrays number of rays within the compact spherical 01406 representation 01407 nnz number of voxels within the sphere 01408 dm an array of size 9 storing transformation 01409 associated with the projection direction 01410 origin coordinates of the center of the volume 01411 ri radius of the sphere 01412 ptrs the beginning address of each ray 01413 cord the coordinates of the first point in each ray 01414 x 3d input volume 01415 y 2d output image 01416 */ 01417 01418 int iqx, iqy, i, j, xc, yc, zc; 01419 float ct, dipx, dipy, dipx1m, dipy1m, xb, yb, dm1, dm4; 01420 int status = 0; 01421 01422 int xcent = origin[0]; 01423 int ycent = origin[1]; 01424 int zcent = origin[2]; 01425 01426 int nx = volsize[0]; 01427 01428 dm1 = dm(1); 01429 dm4 = dm(4); 01430 01431 if ( nx > 2*ri ) { 01432 for (i = 1; i <= nrays; i++) { 01433 zc = cord(1,i)-zcent; 01434 yc = cord(2,i)-ycent; 01435 xc = cord(3,i)-xcent; 01436 01437 xb = zc*dm(1)+yc*dm(2)+xc*dm(3) + xcent; 01438 yb = zc*dm(4)+yc*dm(5)+xc*dm(6) + ycent; 01439 01440 for (j = ptrs(i); j< ptrs(i+1); j++) { 01441 iqx = ifix(xb); 01442 iqy = ifix(yb); 01443 01444 ct = x(j); 01445 dipx = xb - (float)(iqx); 01446 dipy = (yb - (float)(iqy)) * ct; 01447 01448 dipy1m = ct - dipy; 01449 dipx1m = 1.0f - dipx; 01450 01451 y(iqx ,iqy) = y(iqx ,iqy) + dipx1m*dipy1m; 01452 y(iqx+1,iqy) = y(iqx+1,iqy) + dipx*dipy1m; 01453 y(iqx+1,iqy+1) = y(iqx+1,iqy+1) + dipx*dipy; 01454 y(iqx ,iqy+1) = y(iqx ,iqy+1) + dipx1m*dipy; 01455 01456 xb += dm1; 01457 yb += dm4; 01458 } 01459 } 01460 } 01461 else { 01462 fprintf(stderr, " nx must be greater than 2*ri\n"); 01463 exit(1); 01464 } 01465 return status; 01466 }
string EMAN::ChaoProjector::get_desc | ( | ) | const [inline, virtual] |
string EMAN::ChaoProjector::get_name | ( | ) | const [inline, virtual] |
Get the projector's name.
Each projector is indentified by unique name.
Implements EMAN::Projector.
Definition at line 344 of file projector.h.
References NAME.
00345 { 00346 return NAME; 00347 }
TypeDict EMAN::ChaoProjector::get_param_types | ( | ) | const [inline, virtual] |
Get processor parameter information in a dictionary.
Each parameter has one record in the dictionary. Each record contains its name, data-type, and description.
Reimplemented from EMAN::Projector.
Definition at line 359 of file projector.h.
References EMAN::EMObject::FLOAT, EMAN::EMObject::FLOATARRAY, EMAN::EMObject::INT, EMAN::TypeDict::put(), and EMAN::EMObject::TRANSFORM.
00360 { 00361 TypeDict d; 00362 d.put("transform", EMObject::TRANSFORM); 00363 d.put("origin_x", EMObject::INT); 00364 d.put("origin_y", EMObject::INT); 00365 d.put("origin_z", EMObject::INT); 00366 d.put("anglelist", EMObject::FLOATARRAY); 00367 d.put("radius", EMObject::FLOAT); 00368 return d; 00369 }
int ChaoProjector::getnnz | ( | Vec3i | volsize, | |
int | ri, | |||
Vec3i | origin, | |||
int * | nray, | |||
int * | nnz | |||
) | const [private] |
Definition at line 1244 of file projector.cpp.
References nnz, nrays, nx, ny, and status.
Referenced by backproject3d(), and project3d().
01246 : count the number of voxels within a sphere centered 01247 at origin and with a radius ri. 01248 01249 input: 01250 volsize contains the size information (nx,ny,nz) about the volume 01251 ri radius of the object embedded in the cube. 01252 origin coordinates for the center of the volume 01253 01254 output: 01255 nnz total number of voxels within the sphere (of radius ri) 01256 nrays number of rays in z-direction. 01257 */ 01258 { 01259 int ix, iy, iz, rs, r2, xs, ys, zs, xx, yy, zz; 01260 int ftm=0, status = 0; 01261 01262 r2 = ri*ri; 01263 *nnz = 0; 01264 *nrays = 0; 01265 int nx = (int)volsize[0]; 01266 int ny = (int)volsize[1]; 01267 int nz = (int)volsize[2]; 01268 01269 int xcent = (int)origin[0]; 01270 int ycent = (int)origin[1]; 01271 int zcent = (int)origin[2]; 01272 01273 // need to add some error checking 01274 for (ix = 1; ix <=nx; ix++) { 01275 xs = ix-xcent; 01276 xx = xs*xs; 01277 for (iy = 1; iy <= ny; iy++) { 01278 ys = iy-ycent; 01279 yy = ys*ys; 01280 ftm = 1; 01281 for (iz = 1; iz <= nz; iz++) { 01282 zs = iz-zcent; 01283 zz = zs*zs; 01284 rs = xx + yy + zz; 01285 if (rs <= r2) { 01286 (*nnz)++; 01287 if (ftm) { 01288 (*nrays)++; 01289 ftm = 0; 01290 } 01291 } 01292 } 01293 } // end for iy 01294 } // end for ix 01295 return status; 01296 }
int ChaoProjector::ifix | ( | float | a | ) | const [private] |
Definition at line 1547 of file projector.cpp.
Referenced by bckpj3(), and fwdpj3().
01548 { 01549 int ia; 01550 01551 if (a>=0) { 01552 ia = (int)floor(a); 01553 } 01554 else { 01555 ia = (int)ceil(a); 01556 } 01557 return ia; 01558 }
static Projector* EMAN::ChaoProjector::NEW | ( | ) | [inline, static] |
Project an 3D image into a 2D image.
Implements EMAN::Projector.
Definition at line 1602 of file projector.cpp.
References anglelist, cb2sph(), cord, cube, dm, EMDeleteArray(), fwdpj3(), EMAN::EMData::get_data(), EMAN::Util::get_min(), EMAN::Transform::get_rotation(), EMAN::EMData::get_xsize(), EMAN::EMData::get_ysize(), EMAN::EMData::get_zsize(), getnnz(), EMAN::Dict::has_key(), images, LOGERR, nnz, nrays, NullPointerException, EMAN::Projector::params, phi, ptrs, EMAN::EMData::set_attr(), EMAN::EMData::set_complex(), EMAN::EMData::set_ri(), EMAN::EMData::set_size(), setdm(), EMAN::Dict::size(), sphere, status, theta, and EMAN::EMData::update().
01603 { 01604 01605 int nrays, nnz, status, j; 01606 float *dm; 01607 int *ptrs, *cord; 01608 float *sphere, *images; 01609 01610 int nxvol = vol->get_xsize(); 01611 int nyvol = vol->get_ysize(); 01612 int nzvol = vol->get_zsize(); 01613 Vec3i volsize(nxvol,nyvol,nzvol); 01614 01615 int dim = Util::get_min(nxvol,nyvol,nzvol); 01616 if (nzvol == 1) { 01617 LOGERR("The ChaoProjector needs a volume!"); 01618 return 0; 01619 } 01620 Vec3i origin(0,0,0); 01621 // If a sensible origin isn't passed in, choose the middle of 01622 // the cube. 01623 if (params.has_key("origin_x")) {origin[0] = params["origin_x"];} 01624 else {origin[0] = nxvol/2+1;} 01625 if (params.has_key("origin_y")) {origin[1] = params["origin_y"];} 01626 else {origin[1] = nyvol/2+1;} 01627 if (params.has_key("origin_z")) {origin[2] = params["origin_z"];} 01628 else {origin[2] = nzvol/2+1;} 01629 01630 int ri; 01631 if (params.has_key("radius")) {ri = params["radius"];} 01632 else {ri = dim/2 - 1;} 01633 01634 // retrieve the voxel values 01635 float *cube = vol->get_data(); 01636 01637 // count the number of voxels within a sphere centered at icent, 01638 // with radius ri 01639 status = getnnz(volsize, ri, origin, &nrays, &nnz); 01640 // need to check status... 01641 01642 // convert from cube to sphere 01643 sphere = new float[nnz]; 01644 ptrs = new int[nrays+1]; 01645 cord = new int[3*nrays]; 01646 if (sphere == NULL || ptrs == NULL || cord == NULL) { 01647 fprintf(stderr,"ChaoProjector::project3d, failed to allocate!\n"); 01648 exit(1); 01649 } 01650 for (int i = 0; i<nnz; i++) sphere[i] = 0.0; 01651 for (int i = 0; i<nrays+1; i++) ptrs[i] = 0; 01652 for (int i = 0; i<3*nrays; i++) cord[i] = 0; 01653 01654 status = cb2sph(cube, volsize, ri, origin, nnz, ptrs, cord, sphere); 01655 // check status 01656 01657 int nangles = 0; 01658 vector<float> anglelist; 01659 string angletype = "SPIDER"; 01660 // Do we have a list of angles? 01661 if (params.has_key("anglelist")) { 01662 anglelist = params["anglelist"]; 01663 nangles = anglelist.size() / 3; 01664 } else { 01665 Transform* t3d = params["transform"]; 01666 if ( t3d == NULL ) throw NullPointerException("The transform object (required for projection), was not specified"); 01667 // This part was modified by David Woolford - 01668 // Before this the code worked only for SPIDER and EMAN angles, 01669 // but the framework of the Transform3D allows for a generic implementation 01670 // as specified here. 01671 Dict p = t3d->get_rotation("spider"); 01672 if(t3d) {delete t3d; t3d=0;} 01673 01674 float phi = p["phi"]; 01675 float theta = p["theta"]; 01676 float psi = p["psi"]; 01677 anglelist.push_back(phi); 01678 anglelist.push_back(theta); 01679 anglelist.push_back(psi); 01680 nangles = 1; 01681 } 01682 // End David Woolford modifications 01683 01684 dm = new float[nangles*9]; 01685 setdm(anglelist, angletype, dm); 01686 01687 // return images 01688 EMData *ret = new EMData(); 01689 ret->set_size(nxvol, nyvol, nangles); 01690 ret->set_complex(false); 01691 ret->set_ri(true); 01692 01693 images = ret->get_data(); 01694 01695 for (j = 1; j <= nangles; j++) { 01696 status = fwdpj3(volsize, nrays, nnz , &dm(1,j), origin, ri, 01697 ptrs , cord, sphere, &images(1,1,j)); 01698 // check status? 01699 } 01700 01701 // deallocate all temporary work space 01702 EMDeleteArray(dm); 01703 EMDeleteArray(ptrs); 01704 EMDeleteArray(cord); 01705 EMDeleteArray(sphere); 01706 01707 if (!params.has_key("anglelist")) { 01708 Transform* t3d = params["transform"]; 01709 ret->set_attr("xform.projection",t3d); 01710 if(t3d) {delete t3d; t3d=0;} 01711 } 01712 ret->update(); 01713 return ret; 01714 }
void ChaoProjector::setdm | ( | vector< float > | anglelist, | |
string const | angletype, | |||
float * | dm | |||
) | const [private] |
Definition at line 1564 of file projector.cpp.
References anglelist, dgr_to_rad, dm, phi, and theta.
Referenced by backproject3d(), and project3d().
01565 { // convert Euler angles to transformations, dm is an 9 by nangles array 01566 01567 float psi, theta, phi; 01568 double cthe, sthe, cpsi, spsi, cphi, sphi; 01569 int j; 01570 01571 int nangles = anglelist.size() / 3; 01572 01573 // now convert all angles 01574 for (j = 1; j <= nangles; j++) { 01575 phi = static_cast<float>(anglelist(1,j)*dgr_to_rad); 01576 theta = static_cast<float>(anglelist(2,j)*dgr_to_rad); 01577 psi = static_cast<float>(anglelist(3,j)*dgr_to_rad); 01578 01579 // cout << phi << " " << theta << " " << psi << endl; 01580 cthe = cos(theta); 01581 sthe = sin(theta); 01582 cpsi = cos(psi); 01583 spsi = sin(psi); 01584 cphi = cos(phi); 01585 sphi = sin(phi); 01586 01587 dm(1,j)=static_cast<float>(cphi*cthe*cpsi-sphi*spsi); 01588 dm(2,j)=static_cast<float>(sphi*cthe*cpsi+cphi*spsi); 01589 dm(3,j)=static_cast<float>(-sthe*cpsi); 01590 dm(4,j)=static_cast<float>(-cphi*cthe*spsi-sphi*cpsi); 01591 dm(5,j)=static_cast<float>(-sphi*cthe*spsi+cphi*cpsi); 01592 dm(6,j)=static_cast<float>(sthe*spsi); 01593 dm(7,j)=static_cast<float>(sthe*cphi); 01594 dm(8,j)=static_cast<float>(sthe*sphi); 01595 dm(9,j)=static_cast<float>(cthe); 01596 } 01597 }
int ChaoProjector::sph2cb | ( | float * | sphere, | |
Vec3i | volsize, | |||
int | nray, | |||
int | ri, | |||
int | nnz0, | |||
int * | ptrs, | |||
int * | cord, | |||
float * | cube | |||
) | const [private] |
Definition at line 1363 of file projector.cpp.
References cord, cube, nnz, nx, ny, ptrs, sphere, and status.
Referenced by backproject3d().
01365 { 01366 int status=0; 01367 int r2, i, j, ix, iy, iz, nnz; 01368 01369 int nx = (int)volsize[0]; 01370 int ny = (int)volsize[1]; 01371 // int nz = (int)volsize[2]; 01372 01373 r2 = ri*ri; 01374 nnz = 0; 01375 ptrs(1) = 1; 01376 01377 // no need to initialize 01378 // for (i = 0; i<nx*ny*nz; i++) cube[i]=0.0; 01379 01380 nnz = 0; 01381 for (j = 1; j <= nrays; j++) { 01382 iz = cord(1,j); 01383 iy = cord(2,j); 01384 ix = cord(3,j); 01385 for (i = ptrs(j); i<=ptrs(j+1)-1; i++, iz++) { 01386 nnz++; 01387 cube(iz,iy,ix) = sphere(nnz); 01388 } 01389 } 01390 if (nnz != nnz0) status = -1; 01391 return status; 01392 }
const string ChaoProjector::NAME = "chao" [static] |