fundamentals.cpp

Go to the documentation of this file.
00001 /*
00002  * Author: Pawel A.Penczek, 09/09/2006 (Pawel.A.Penczek@uth.tmc.edu)
00003  * Copyright (c) 2000-2006 The University of Texas - Houston Medical School
00004  *
00005  * This software is issued under a joint BSD/GNU license. You may use the
00006  * source code in this file under either license. However, note that the
00007  * complete EMAN2 and SPARX software packages have some GPL dependencies,
00008  * so you are responsible for compliance with the licenses of these packages
00009  * if you opt to use BSD licensing. The warranty disclaimer below holds
00010  * in either instance.
00011  *
00012  * This complete copyright notice must be included in any revised version of the
00013  * source code. Additional authorship citations may be added, but existing
00014  * author citations must be preserved.
00015  *
00016  * This program is free software; you can redistribute it and/or modify
00017  * it under the terms of the GNU General Public License as published by
00018  * the Free Software Foundation; either version 2 of the License, or
00019  * (at your option) any later version.
00020  *
00021  * This program is distributed in the hope that it will be useful,
00022  * but WITHOUT ANY WARRANTY; without even the implied warranty of
00023  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
00024  * GNU General Public License for more details.
00025  *
00026  * You should have received a copy of the GNU General Public License
00027  * along with this program; if not, write to the Free Software
00028  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
00029  *
00030  */
00031 
00032 #include "emdata.h"
00033 
00034 using namespace EMAN;
00035 using std::complex;
00036 
00037 namespace EMAN {
00038 
00039 EMData* 
00040 periodogram(EMData* f) {
00041         // These are actual dimensions
00042         int nx  = f->get_xsize();
00043         int ny  = f->get_ysize();
00044         int nz  = f->get_zsize();
00045                 // We manifestly assume no zero-padding here, just the 
00046                 // necessary extension along x for the fft
00047 
00048         if (f->is_complex()) nx = (nx - 2 + f->is_fftodd()); // nx is the real-space size of the input image
00049         int lsd2 = (nx + 2 - nx%2) / 2; // Extended x-dimension of the complex image
00050 
00051 //  Process f if real
00052         EMData* fp = NULL;
00053         if(f->is_complex()) fp = f->copy(); // we need to make a full copy so that we don't damage the original
00054         else {
00055                 fp = f->norm_pad(false, 1); // Extend and do the FFT if f is real
00056                 fp->do_fft_inplace();
00057         }
00058         fp->set_array_offsets(1,1,1);
00059 
00060         //  Periodogram: fp:=|fp|**2
00061         for (int iz = 1; iz <= nz; iz++) {
00062                 for (int iy = 1; iy <= ny; iy++) {
00063                         for (int ix = 1; ix <= lsd2; ix++) {
00064                                 float fpr = real(fp->cmplx(ix,iy,iz));
00065                                 float fpi = imag(fp->cmplx(ix,iy,iz));
00066                                 fp->cmplx(ix,iy,iz) = fpr*fpr + fpi*fpi;
00067                         }
00068                 }
00069         }
00070         //  Create power as a 3D array (-n/2:n/2+n%2-1)
00071         int nyt, nzt;
00072         int nx2 = nx/2;
00073         int ny2 = ny/2; if(ny2 == 0) nyt =0; else nyt=ny;
00074         int nz2 = nz/2; if(nz2 == 0) nzt =0; else nzt=nz;
00075         int nx2p = nx2+nx%2;
00076         int ny2p = ny2+ny%2;
00077         int nz2p = nz2+nz%2;
00078         EMData& power = *(new EMData()); // output image
00079         power.set_size(nx, ny, nz);
00080         power.set_array_offsets(-nx2,-ny2,-nz2);
00081 //If instead of preservation of the norm one would prefer to have peak of a PW of a single sine wave equal one
00082 //                             multiply power by the scale below, or the other way around.
00083         float scale = 4.0f/float (nx*nx)/float (ny*ny)/float (nz*nz);
00084         for (int iz = 1; iz <= nz; iz++) {
00085                 int jz=iz-1; 
00086                 if(jz>=nz2p) jz=jz-nzt;
00087                 for (int iy = 1; iy <= ny; iy++) {
00088                         int jy=iy-1; 
00089                         if(jy>=ny2p) jy=jy-nyt;
00090                         for (int ix = 1; ix <= lsd2; ix++) {
00091                                 int jx=ix-1; 
00092                                 if(jx>=nx2p) jx=jx-nx;
00093                                 power(jx,jy,jz) = real(fp->cmplx(ix,iy,iz)) * scale;
00094                         }
00095                 }
00096         }
00097 
00098 //  Create the Friedel related half
00099         int  nzb, nze, nyb, nye, nxb, nxe;
00100         nxb =-nx2+(nx+1)%2;
00101         nxe = nx2-(nx+1)%2;
00102         if(ny2 == 0) {nyb =0; nye = 0;} else {nyb =-ny2+(ny+1)%2; nye = ny2-(ny+1)%2;}
00103         if(nz2 == 0) {nzb =0; nze = 0;} else {nzb =-nz2+(nz+1)%2; nze = nz2-(nz+1)%2;}
00104         for (int iz = nzb; iz <= nze; iz++) {
00105                 for (int iy = nyb; iy <= nye; iy++) {
00106                         for (int ix = 1; ix <= nxe; ix++) { // Note this loop begins with 1 - FFT should create correct Friedel related 0 plane
00107                                 power(-ix,-iy,-iz) = power(ix,iy,iz);
00108                         }
00109                 }
00110         }
00111         if(ny2 != 0)  {
00112                 if(nz2 != 0)  {
00113                         if(nz%2 == 0) {  //if nz even, fix the first slice
00114                                 for (int iy = nyb; iy <= nye; iy++) {
00115                                         for (int ix = nxb; ix <= -1; ix++) { 
00116                                                 power(ix,iy,-nz2) = power(-ix,-iy,-nz2);
00117                                         }
00118                                 }
00119                                 if(ny%2 == 0) {  //if ny even, fix the first line
00120                                         for (int ix = nxb; ix <= -1; ix++) { 
00121                                                 power(ix,-ny2,-nz2) = power(-ix,-ny2,-nz2);
00122                                         }
00123                                 }
00124                         }
00125                 }
00126                 if(ny%2 == 0) {  //if ny even, fix the first column
00127                         for (int iz = nzb; iz <= nze; iz++) {
00128                                 for (int ix = nxb; ix <= -1; ix++) {
00129                                         power(ix,-ny2,-iz) = power(-ix,-ny2,iz);
00130                                 }
00131                         }
00132                 }
00133                 
00134         }
00135                 
00136         if( fp ) {
00137                 delete fp; // avoid a memory leak!
00138                 fp = 0;
00139         }
00140         //power[0][0][0]=power[1][0][0];  //Steve requested the original origin.
00141 
00142         power.update();
00143         power.set_array_offsets(0,0,0);
00144         return &power;
00145 //OVER AND OUT
00146 }
00147 
00148 
00149 /*
00150 fourierproduct
00151 Purpose: Calculate various fundamental image processing operations 
00152          based on Fourier products for 1-2-3D images. 
00153 Method: Calculate FFT (if needed),
00154         Real f input -> padded fp workspace
00155         Complex f input -> set fp to f
00156         Real g input -> padded gp workspace
00157         Complex g input -> set gp to g
00158         No g input -> set gp tp fp
00159          -> real output -> if f real delete workspace fp and if g real delete workspace gp.
00160 Input: f real or complex 1-2-3D image, g - real or complex 1-2-3D image.
00161 Output: 1-2-3D real image with the result
00162         (correlation F*conjg(G), convolution F*G, 
00163                 autocorrelation |F|^2, self-correlation |F|)
00164         and with the origin at the image center int[n/2]+1.
00165 */  
00166         EMData* 
00167         fourierproduct(EMData* f, EMData* g, fp_flag flag, fp_type ptype, bool center) {
00168                 int normfact;
00169                 //std::complex<float> phase_mult;
00170                 // Not only does the value of "flag" determine how we handle
00171                 // periodicity, but it also determines whether or not we should
00172                 // normalize the results.  Here's some convenience bools:
00173                 bool donorm = (0 == flag%2) ? true : false;
00174                 // the 2x padding is hardcoded for now
00175                 int  npad  = (flag >= 3) ? 2 : 1;  // amount of padding used
00176                 // g may be NULL.  If so, have g point to the same object as f.  In that
00177                 // case we need to be careful later on not to try to delete g's workspace
00178                 // as well as f's workspace, since they will be the same.
00179                 bool  gexists = true;
00180                 if (!g) { g = f; gexists = false; }
00181                 if ( f->is_complex() || g->is_complex() ) {
00182                         // Fourier input only allowed for circulant
00183                         if (CIRCULANT != flag) {
00184                                 LOGERR("Cannot perform normalization or padding on Fourier type.");
00185                                 throw InvalidValueException(flag, "Cannot perform normalization or padding on Fourier type.");
00186                         }
00187                 }
00188                 // These are actual dimensions of f (and real-space sizes for ny and nz)
00189                 int nx  = f->get_xsize();
00190                 int ny  = f->get_ysize();
00191                 int nz  = f->get_zsize();
00192                 // We manifestly assume no zero-padding here, just the 
00193                 // necessary extension along x for the fft
00194                 if (!f->is_real()) nx = (nx - 2 + (f->is_fftodd() ? 1 : 0)); 
00195 
00196                 // these are padded dimensions
00197                 const int nxp = npad*nx;
00198                 const int nyp = (ny > 1) ? npad*ny : 1; // don't pad y for 1-d image
00199                 const int nzp = (nz > 1) ? npad*nz : 1; // don't pad z for 2-d image
00200 
00201                 // now one half of the padded, fft-extended size along x
00202                 const int lsd2 = (nxp + 2 - nxp%2) / 2; 
00203                 // The [padded] fft-extended fourier version of f is fp.
00204 
00205                 EMData* fp = NULL;
00206                 if (f->is_complex()) { 
00207                         // If f is already a fourier object then fp is a copy of f.
00208                         // (The fp workspace is modified, so we copy f to keep f pristine.)
00209                         fp=f->copy();
00210                 } else {
00211                         //  [normalize] [pad] compute fft
00212                         fp = f->norm_pad(donorm, npad);
00213                         fp->do_fft_inplace();
00214                 }
00215                 // The [padded] fft-extended version of g is gp.
00216                 EMData* gp = NULL;
00217                 if(f==g) {
00218                         // g is an alias for f, so gp should be an alias for fp
00219                         gp=fp;
00220                 } else if (g->is_complex()) {
00221                         // g is already a Fourier object, so gp is just an alias for g
00222                         // (The gp workspace is not modified, so we don't need a copy.)
00223                         gp = g;
00224                 } else {
00225                         // normal case: g is real and different from f, so compute gp
00226                         gp = g->norm_pad(donorm, npad);
00227                         gp->do_fft_inplace();
00228                 }
00229                 // Get complex matrix views of fp and gp; matrices start from 1 (not 0)
00230                 fp->set_array_offsets(1,1,1);
00231                 gp->set_array_offsets(1,1,1);
00232                 
00233                 // If the center flag is true, put the center of the correlation in the middle
00234                 // If it is false, put it in (0,0), this approach saves time, but it is diffcult to manage the result
00235                 if (center) {
00236                         //  Multiply two functions (the real work of this routine)
00237                         int itmp = nx/2;
00238                         //float sx  = float(-twopi*float(itmp)/float(nxp));
00239                         float sxn = 2*float(itmp)/float(nxp);
00240                         float sx = -M_PI*sxn;
00241                         itmp = ny/2;
00242                         //float sy  = float(-twopi*float(itmp)/float(nyp));
00243                         float syn = 2*float(itmp)/float(nyp);
00244                         float sy = -M_PI*syn;
00245                         itmp = nz/2;
00246                         //float sz  = float(-twopi*float(itmp)/float(nzp));
00247                         float szn = 2*float(itmp)/float(nzp);
00248                         float sz = -M_PI*szn;
00249                         if ( nx%2==0 && (ny%2==0 || ny==1 ) && (nz%2==0 || nz==1 ) ) {
00250                                 switch (ptype) {
00251                                         case AUTOCORRELATION:
00252                                         // fpmat := |fpmat|^2
00253                                         // Note nxp are padded dimensions
00254                                                 for (int iz = 1; iz <= nzp; iz++) {
00255                                                         for (int iy = 1; iy <= nyp; iy++) {
00256                                                                 for (int ix = 1; ix <= lsd2; ix++) {
00257                                                                         float fpr = real(fp->cmplx(ix,iy,iz));
00258                                                                         float fpi = imag(fp->cmplx(ix,iy,iz));
00259                                                                         fp->cmplx(ix,iy,iz) = complex<float>(fpr*fpr+fpi*fpi, 0.0f);
00260                                                                 }
00261                                                         }
00262                                                 }
00263                                                 break;
00264                                         case SELF_CORRELATION:
00265                                         // fpmat:=|fpmat|
00266                                         // Note nxp are padded dimensions
00267                                                 for (int iz = 1; iz <= nzp; iz++) {
00268                                                         for (int iy = 1; iy <= nyp; iy++) {
00269                                                                 for (int ix = 1; ix <= lsd2; ix++) {
00270                                                                         fp->cmplx(ix,iy,iz) = complex<float>(abs(fp->cmplx(ix,iy,iz)), 0.0f);
00271                                                                 }
00272                                                         }
00273                                                 }
00274                                                 break;
00275                                         case CORRELATION:
00276                                         // fpmat:=fpmat*conjg(gpmat)
00277                                         // Note nxp are padded dimensions
00278                                                 for (int iz = 1; iz <= nzp; iz++) {
00279                                                         for (int iy = 1; iy <= nyp; iy++) {
00280                                                                 for (int ix = 1; ix <= lsd2; ix++) {
00281                                                                         fp->cmplx(ix,iy,iz) *= conj(gp->cmplx(ix,iy,iz));
00282                                                                 }
00283                                                         }
00284                                                 }
00285                                         break;
00286                                         case CONVOLUTION:
00287                                         // fpmat:=fpmat*gpmat
00288                                         // Note nxp are padded dimensions
00289                                                 for (int iz = 1; iz <= nzp; iz++) {
00290                                                         for (int iy = 1; iy <= nyp; iy++) {
00291                                                                 for (int ix = 1; ix <= lsd2; ix++) {
00292                                                                         fp->cmplx(ix,iy,iz) *= gp->cmplx(ix,iy,iz);
00293                                                                 }
00294                                                         }
00295                                                 }
00296                                                 break;
00297                                         default:
00298                                                 LOGERR("Illegal option in Fourier Product");
00299                                                 throw InvalidValueException(ptype, "Illegal option in Fourier Product");
00300                                 }                                       
00301                                 for (int iz = 1; iz <= nzp; iz++) {
00302                                         for (int iy = 1; iy <= nyp; iy++) {
00303                                                 for (int ix = (iz+iy+1)%2+1; ix <= lsd2; ix+=2) {
00304                                                         fp->cmplx(ix,iy,iz) = -fp->cmplx(ix,iy,iz);
00305                                                 }
00306                                         }
00307                                 }
00308                         } else {                        
00309                                 switch (ptype) {
00310                                         case AUTOCORRELATION:
00311                                         // fpmat := |fpmat|^2
00312                                         // Note nxp are padded dimensions
00313                                                 for (int iz = 1; iz <= nzp; iz++) {
00314                                                 int jz=iz-1; if(jz>nzp/2) jz=jz-nzp; float argz=sz*jz;
00315                                                         for (int iy = 1; iy <= nyp; iy++) {
00316                                                         int jy=iy-1; if(jy>nyp/2) jy=jy-nyp; float argy=sy*jy+argz;
00317                                                                 for (int ix = 1; ix <= lsd2; ix++) {
00318                                                                         int jx=ix-1; float arg=sx*jx+argy;
00319                                                                         float fpr = real(fp->cmplx(ix,iy,iz));
00320                                                                         float fpi = imag(fp->cmplx(ix,iy,iz));
00321                                                                         fp->cmplx(ix,iy,iz)= (fpr*fpr + fpi*fpi) *std::complex<float>(cos(arg),sin(arg));
00322                                                                 }
00323                                                         }
00324                                                 }
00325                                                 break;
00326                                         case SELF_CORRELATION:
00327                                         // fpmat:=|fpmat|
00328                                         // Note nxp are padded dimensions
00329                                                 for (int iz = 1; iz <= nzp; iz++) {
00330                                                 int jz=iz-1; if(jz>nzp/2) jz=jz-nzp; float argz=sz*jz;
00331                                                         for (int iy = 1; iy <= nyp; iy++) {
00332                                                         int jy=iy-1; if(jy>nyp/2) jy=jy-nyp; float argy=sy*jy+argz;
00333                                                                 for (int ix = 1; ix <= lsd2; ix++) {
00334                                                                         int jx=ix-1; float arg=sx*jx+argy;
00335                                                                         fp->cmplx(ix,iy,iz) = abs(fp->cmplx(ix,iy,iz)) *std::complex<float>(cos(arg),sin(arg));
00336                                                                 }
00337                                                         }
00338                                                 }
00339                                                 break;
00340                                         case CORRELATION:
00341                                         // fpmat:=fpmat*conjg(gpmat)
00342                                         // Note nxp are padded dimensions
00343                                                 for (int iz = 1; iz <= nzp; iz++) {
00344                                                 int jz=iz-1; if(jz>nzp/2) jz=jz-nzp; float argz=sz*jz;
00345                                                         for (int iy = 1; iy <= nyp; iy++) {
00346                                                         int jy=iy-1; if(jy>nyp/2) jy=jy-nyp; float argy=sy*jy+argz;
00347                                                                 for (int ix = 1; ix <= lsd2; ix++) {
00348                                                                         int jx=ix-1; float arg=sx*jx+argy;
00349                                                                         fp->cmplx(ix,iy,iz) *= conj(gp->cmplx(ix,iy,iz)) *std::complex<float>(cos(arg),sin(arg));
00350                                                                 }
00351                                                         }
00352                                                 }
00353                                         break;
00354                                         case CONVOLUTION:
00355                                         // fpmat:=fpmat*gpmat
00356                                         // Note nxp are padded dimensions
00357                                                 if(npad == 1) {
00358                                                         sx -= 4*(nx%2)/float(nx);
00359                                                         sy -= 4*(ny%2)/float(ny);
00360                                                         sz -= 4*(nz%2)/float(nz);
00361                                                 }
00362                                                 for (int iz = 1; iz <= nzp; iz++) {
00363                                                         int jz=iz-1; if(jz>nzp/2) jz=jz-nzp; float argz=sz*jz;
00364                                                         for (int iy = 1; iy <= nyp; iy++) {
00365                                                                 int jy=iy-1; if(jy>nyp/2) jy=jy-nyp; float argy=sy*jy+argz;
00366                                                                 for (int ix = 1; ix <= lsd2; ix++) {
00367                                                                         int jx=ix-1; float arg=sx*jx+argy;
00368                                                                         fp->cmplx(ix,iy,iz) *= gp->cmplx(ix,iy,iz) *std::complex<float>(cos(arg),sin(arg));
00369                                                                 }
00370                                                         }
00371                                                 }
00372                                                 break;
00373                                         default:
00374                                                 LOGERR("Illegal option in Fourier Product");
00375                                                 throw InvalidValueException(ptype, "Illegal option in Fourier Product");
00376                                 }
00377                         }
00378                 } else {
00379                         // If the center flag is false, then just do basic multiplication
00380                         switch (ptype) {
00381                                 case AUTOCORRELATION:
00382                                         for (int iz = 1; iz <= nzp; iz++) {
00383                                                 for (int iy = 1; iy <= nyp; iy++) {
00384                                                         for (int ix = 1; ix <= lsd2; ix++) {
00385                                                                 float fpr = real(fp->cmplx(ix,iy,iz));
00386                                                                 float fpi = imag(fp->cmplx(ix,iy,iz));
00387                                                                 fp->cmplx(ix,iy,iz) = complex<float>(fpr*fpr+fpi*fpi, 0.0f);
00388                                                         }
00389                                                 }
00390                                         }
00391                                         break;
00392                                 case SELF_CORRELATION:
00393                                         for (int iz = 1; iz <= nzp; iz++) {
00394                                                 for (int iy = 1; iy <= nyp; iy++) {
00395                                                         for (int ix = 1; ix <= lsd2; ix++) {
00396                                                                 fp->cmplx(ix,iy,iz) = complex<float>(abs(fp->cmplx(ix,iy,iz)), 0.0f);
00397                                                         }
00398                                                 }
00399                                         }
00400                                         break;
00401                                 case CORRELATION:
00402                                         //phase_mult = 1;
00403                                         for (int iz = 1; iz <= nzp; iz++) {
00404                                                 for (int iy = 1; iy <= nyp; iy++) {
00405                                                         for (int ix = 1; ix <= lsd2; ix++) {
00406                                                                 fp->cmplx(ix,iy,iz)*= conj(gp->cmplx(ix,iy,iz));
00407                                                         }
00408                                                 }
00409                                         }
00410                                         break;
00411                                 case CONVOLUTION:
00412                                         if(npad == 1) {
00413                                                 float sx = -M_PI*2*(nx%2)/float(nx);
00414                                                 float sy = -M_PI*2*(ny%2)/float(ny);
00415                                                 float sz = -M_PI*2*(nz%2)/float(nz);
00416                                                 for (int iz = 1; iz <= nzp; iz++) {
00417                                                         int jz=iz-1; if(jz>nzp/2) jz=jz-nzp; float argz=sz*jz;
00418                                                         for (int iy = 1; iy <= nyp; iy++) {
00419                                                                 int jy=iy-1; if(jy>nyp/2) jy=jy-nyp; float argy=sy*jy+argz;
00420                                                                 for (int ix = 1; ix <= lsd2; ix++) {
00421                                                                         int jx=ix-1; float arg=sx*jx+argy;
00422                                                                         fp->cmplx(ix,iy,iz) *= gp->cmplx(ix,iy,iz) *std::complex<float>(cos(arg),sin(arg));
00423                                                                 }
00424                                                         }
00425                                                 }
00426                                         } else {
00427                                                 for (int iz = 1; iz <= nzp; iz++) {
00428                                                         for (int iy = 1; iy <= nyp; iy++) {
00429                                                                 for (int ix = 1; ix <= lsd2; ix++) {
00430                                                                         fp->cmplx(ix,iy,iz)*= gp->cmplx(ix,iy,iz);
00431                                                                 }
00432                                                         }
00433                                                 }
00434                                         }
00435                                         break;
00436                                 default:
00437                                         LOGERR("Illegal option in Fourier Product");
00438                                         throw InvalidValueException(ptype, "Illegal option in Fourier Product");
00439                         }
00440                 }
00441                 // Now done w/ gp, so let's get rid of it (if it's not an alias of fp or simply g was complex on input);
00442                 if (gexists && (f != g) && (!g->is_complex())) {
00443                         if( gp ) {
00444                                 delete gp;
00445                                 gp = 0;
00446                         }
00447                 }
00448                 // back transform
00449                 fp->do_ift_inplace();
00450                 if(center && npad ==2)  fp->depad();
00451                 else                    fp->depad_corner();
00452 
00453                 //vector<int> saved_offsets = fp->get_array_offsets();  I do not know what the meaning of it was, did not work anyway PAP
00454                 fp->set_array_offsets(1,1,1);
00455 
00456                 normfact = (nxp/nx)*(nyp/ny)*(nzp/nz);  // Normalization factor for the padded operations
00457                 if(normfact>1) {
00458                         for (int iz = 1; iz <= nz; iz++) for (int iy = 1; iy <= ny; iy++) for (int ix = 1; ix <= nx; ix++) (*fp)(ix,iy,iz) *= normfact;
00459                 }
00460                 // Lag normalization
00461                 if(flag>4)  {
00462                         normfact = nx*ny*nz;  // Normalization factor
00463                         int nxc=nx/2+1, nyc=ny/2+1, nzc=nz/2+1;
00464                         for (int iz = 1; iz <= nz; iz++) {
00465                                 float lagz=float(normfact/(nz-abs(iz-nzc)));
00466                                 for (int iy = 1; iy <= ny; iy++) {
00467                                         float lagyz=lagz/(ny-abs(iy-nyc));
00468                                         for (int ix = 1; ix <= nx; ix++) {
00469                                                 (*fp)(ix,iy,iz) *= lagyz/(nx-abs(ix-nxc));
00470                                         }
00471                                 }
00472                         }       
00473                 }
00474                 //OVER AND OUT
00475                 //fp->set_array_offsets(saved_offsets);  This was strange and did not work, PAP
00476                 fp->set_array_offsets(0,0,0);
00477                 fp->update();
00478                 return fp;
00479         }
00480 }
00481 

Generated on Mon Jul 19 12:40:10 2010 for EMAN2 by  doxygen 1.4.7