EMAN::ChaoProjector Class Reference

Fast real space projection using Bi-Linear interpolation. More...

#include <projector.h>

Inheritance diagram for EMAN::ChaoProjector:

Inheritance graph
[legend]
Collaboration diagram for EMAN::ChaoProjector:

Collaboration graph
[legend]
List of all members.

Public Member Functions

EMDataproject3d (EMData *vol) const
 Project an 3D image into a 2D image.
EMDatabackproject3d (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 ProjectorNEW ()

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

Detailed Description

Fast real space projection using Bi-Linear interpolation.

(C. Yang)

Definition at line 338 of file projector.h.


Member Function Documentation

EMData * ChaoProjector::backproject3d ( EMData imagestack  )  const [virtual]

Back-project a 2D image into a 3D image.

Returns:
A 3D image from the backprojection.

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]

Implements EMAN::Projector.

Definition at line 349 of file projector.h.

00350                 {
00351                         return "Fast real space projection generation with Bi-Linear interpolation.";
00352                 }

string EMAN::ChaoProjector::get_name (  )  const [inline, virtual]

Get the projector's name.

Each projector is indentified by unique name.

Returns:
The projector's 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.

Returns:
A dictionary containing the parameter info.

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]

Definition at line 354 of file projector.h.

00355                 {
00356                         return new ChaoProjector();
00357                 }

EMData * ChaoProjector::project3d ( EMData vol  )  const [virtual]

Project an 3D image into a 2D image.

Returns:
A 2D image from the projection.

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 }


Member Data Documentation

const string ChaoProjector::NAME = "chao" [static]

Definition at line 371 of file projector.h.

Referenced by get_name().


The documentation for this class was generated from the following files:
Generated on Thu Nov 17 12:47:10 2011 for EMAN2 by  doxygen 1.4.7