Xmipp  v3.23.11-Nereus
Public Types | Public Member Functions | Public Attributes | List of all members
ProgForwardArtZernike3D Class Reference

#include <forward_art_zernike3d.h>

Inheritance diagram for ProgForwardArtZernike3D:
Inheritance graph
[legend]
Collaboration diagram for ProgForwardArtZernike3D:
Collaboration graph
[legend]

Public Types

enum  Mode { Mode::Proj, Mode::Vol }
 

Public Member Functions

 ProgForwardArtZernike3D ()
 Empty constructor. More...
 
 ~ProgForwardArtZernike3D ()
 Destructor. More...
 
void readParams ()
 Read argument from command line. More...
 
void show ()
 Show. More...
 
void defineParams ()
 Define parameters. More...
 
void preProcess ()
 
void processImage (const FileName &fnImg, const FileName &fnImgOut, const MDRow &rowIn, MDRow &rowOut)
 
void numCoefficients (int l1, int l2, int &vecSize)
 Length of coefficients vector. More...
 
void fillVectorTerms (int l1, int l2, Matrix1D< int > &vL1, Matrix1D< int > &vN, Matrix1D< int > &vL2, Matrix1D< int > &vM)
 Zernike and SPH coefficients allocation. More...
 
void deformVol (MultidimArray< double > &mP, MultidimArray< double > &mW, const MultidimArray< double > &mV, double rot, double tilt, double psi)
 Deform a volumen using Zernike-Spherical harmonic basis. More...
 
void recoverVol ()
 
virtual void finishProcessing ()
 
double bspline1 (double x)
 
- Public Member Functions inherited from XmippMetadataProgram
MetaDatagetInputMd ()
 
MetaDataVecgetOutputMd ()
 
 XmippMetadataProgram ()
 Empty constructor. More...
 
virtual int tryRead (int argc, const char **argv, bool reportErrors=true)
 
virtual void init ()
 
virtual void setup (MetaData *md, const FileName &o="", const FileName &oroot="", bool applyGeo=false, MDLabel label=MDL_IMAGE)
 
virtual ~XmippMetadataProgram ()
 
void setMode (WriteModeMetaData _mode)
 
void setupRowOut (const FileName &fnImgIn, const MDRow &rowIn, const FileName &fnImgOut, MDRow &rowOut) const
 Prepare rowout. More...
 
virtual void wait ()
 Wait for the distributor to finish. More...
 
virtual void checkPoint ()
 For very long programs, it may be needed to write checkpoints. More...
 
- Public Member Functions inherited from XmippProgram
const char * getParam (const char *param, int arg=0)
 
const char * getParam (const char *param, const char *subparam, int arg=0)
 
int getIntParam (const char *param, int arg=0)
 
int getIntParam (const char *param, const char *subparam, int arg=0)
 
double getDoubleParam (const char *param, int arg=0)
 
double getDoubleParam (const char *param, const char *subparam, int arg=0)
 
float getFloatParam (const char *param, int arg=0)
 
float getFloatParam (const char *param, const char *subparam, int arg=0)
 
void getListParam (const char *param, StringVector &list)
 
int getCountParam (const char *param)
 
bool checkParam (const char *param)
 
bool existsParam (const char *param)
 
void addParamsLine (const String &line)
 
void addParamsLine (const char *line)
 
ParamDefgetParamDef (const char *param) const
 
virtual void quit (int exit_code=0) const
 
virtual int tryRun ()
 
void initProgress (size_t total, size_t stepBin=60)
 
void setProgress (size_t value=0)
 
void endProgress ()
 
void processDefaultComment (const char *param, const char *left)
 
void setDefaultComment (const char *param, const char *comment)
 
void setProgramName (const char *name)
 
void addUsageLine (const char *line, bool verbatim=false)
 
void clearUsage ()
 
void addExampleLine (const char *example, bool verbatim=true)
 
void addSeeAlsoLine (const char *seeAlso)
 
void addKeywords (const char *keywords)
 
const char * name () const
 
virtual void usage (int verb=0) const
 
virtual void usage (const String &param, int verb=2)
 
int version () const
 
virtual void read (int argc, const char **argv, bool reportErrors=true)
 
virtual void read (int argc, char **argv, bool reportErrors=true)
 
void read (const String &argumentsLine)
 
 XmippProgram ()
 
 XmippProgram (int argc, const char **argv)
 
virtual ~XmippProgram ()
 

Public Attributes

FileName fnVolR
 
FileName fnMaskR
 
FileName fnVolO
 
FileName fnOutDir
 Output directory. More...
 
int L1
 
int L2
 
Matrix1D< int > vL1
 
Matrix1D< int > vN
 
Matrix1D< int > vL2
 
Matrix1D< int > vM
 
double Ts
 
int RmaxDef
 
bool phaseFlipped
 
bool ignoreCTF
 
double lambda
 
int save_iter
 
bool useCTF
 
bool useZernike
 
int flagEnabled
 
bool resume
 
int niter
 
int sort_last_N
 
MultidimArray< int > mask2D
 
size_t Xdim
 
Image< double > V
 
Image< double > Vrefined
 
Image< double > Vout
 
Image< double > Ifilteredp
 
Image< double > I
 
MultidimArray< int > Vmask
 
Image< double > P
 
Image< double > W
 
Image< double > Idiff
 
Matrix2D< double > A
 
double rot
 
double tilt
 
double psi
 
double shiftX
 
double shiftY
 
bool flip
 
bool hasCTF
 
double defocusU
 
double defocusV
 
double defocusAngle
 
CTFDescription ctf
 
FourierFilter FilterCTF
 
int vecSize
 
std::vector< double > clnm
 
bool showOptimization
 
MultidimArray< size_t > ordered_list
 
int current_save_iter
 
size_t num_images
 
size_t current_image
 
int current_iter
 
int initX
 
int endX
 
int initY
 
int endY
 
int initZ
 
int endZ
 
int loop_step
 
struct blobtype blob
 
double blob_r
 
double sigma
 
double sigma4
 
Matrix1D< double > gaussianProjectionTable
 
Matrix1D< double > gaussianProjectionTable2
 
FourierFilter filter
 
FourierFilter filter2
 
- Public Attributes inherited from XmippMetadataProgram
FileName fn_in
 Filenames of input and output Metadata. More...
 
FileName fn_out
 
FileName baseName
 
FileName pathBaseName
 
FileName oextBaseName
 
bool apply_geo
 Apply geo. More...
 
size_t ndimOut
 Output dimensions. More...
 
size_t zdimOut
 
size_t ydimOut
 
size_t xdimOut
 
DataType datatypeOut
 
size_t mdInSize
 Number of input elements. More...
 
- Public Attributes inherited from XmippProgram
bool doRun
 
bool runWithoutArgs
 
int verbose
 Verbosity level. More...
 
int debug
 

Additional Inherited Members

- Protected Member Functions inherited from XmippMetadataProgram
virtual void initComments ()
 
virtual void postProcess ()
 
virtual bool getImageToProcess (size_t &objId, size_t &objIndex)
 
void show () const override
 
virtual void startProcessing ()
 
virtual void writeOutput ()
 
virtual void showProgress ()
 
virtual void defineLabelParam ()
 
- Protected Member Functions inherited from XmippProgram
void defineCommons ()
 
- Protected Attributes inherited from XmippMetadataProgram
WriteModeMetaData mode
 Metadata writing mode: OVERWRITE, APPEND. More...
 
FileName oext
 Output extension and root. More...
 
FileName oroot
 
MDLabel image_label
 MDLabel to be used to read/write images, usually will be MDL_IMAGE. More...
 
bool produces_an_output
 Indicate that a unique final output is produced. More...
 
bool produces_a_metadata
 Indicate that the unique final output file is a Metadata. More...
 
bool each_image_produces_an_output
 Indicate that an output is produced for each image in the input. More...
 
bool allow_apply_geo
 
bool decompose_stacks
 Input Metadata will treat a stack file as a set of images instead of a unique file. More...
 
bool delete_output_stack
 Delete previous output stack file prior to process images. More...
 
bool get_image_info
 Get the input image file dimensions to further operations. More...
 
bool save_metadata_stack
 Save the associated output metadata when output file is a stack. More...
 
bool track_origin
 Include the original input image filename in the output stack. More...
 
bool keep_input_columns
 Keep input metadata columns. More...
 
bool remove_disabled
 Remove disabled images from the input selfile. More...
 
bool allow_time_bar
 Show process time bar. More...
 
bool input_is_metadata
 Input is a metadata. More...
 
bool single_image
 Input is a single image. More...
 
bool input_is_stack
 Input is a stack. More...
 
bool output_is_stack
 Output is a stack. More...
 
bool create_empty_stackfile
 
bool delete_mdIn
 
size_t time_bar_step
 Some time bar related counters. More...
 
size_t time_bar_size
 
size_t time_bar_done
 
- Protected Attributes inherited from XmippProgram
int errorCode
 
ProgramDefprogDef
 Program definition and arguments parser. More...
 
std::map< String, CommentListdefaultComments
 
int argc
 Original command line arguments. More...
 
const char ** argv
 

Detailed Description

Predict Continuous Parameters.

Definition at line 40 of file forward_art_zernike3d.h.

Member Enumeration Documentation

◆ Mode

Enumerator
Proj 
Vol 

Definition at line 144 of file forward_art_zernike3d.h.

144 { Proj, Vol };

Constructor & Destructor Documentation

◆ ProgForwardArtZernike3D()

ProgForwardArtZernike3D::ProgForwardArtZernike3D ( )

Empty constructor.

Definition at line 35 of file forward_art_zernike3d.cpp.

36 {
37  resume = false;
38  produces_a_metadata = true;
40  showOptimization = false;
41 }
bool each_image_produces_an_output
Indicate that an output is produced for each image in the input.
bool produces_a_metadata
Indicate that the unique final output file is a Metadata.

◆ ~ProgForwardArtZernike3D()

ProgForwardArtZernike3D::~ProgForwardArtZernike3D ( )
default

Destructor.

Member Function Documentation

◆ bspline1()

double ProgForwardArtZernike3D::bspline1 ( double  x)

Definition at line 759 of file forward_art_zernike3d.cpp.

760 {
761  double m = 1 / sigma;
762  if (0. < x && x < sigma)
763  return m * (sigma - x);
764  else if (-sigma < x && x <= 0.)
765  return m * (sigma + x);
766  else
767  return 0.;
768 }
doublereal * x
int m

◆ defineParams()

void ProgForwardArtZernike3D::defineParams ( )
virtual

Define parameters.

Reimplemented from XmippMetadataProgram.

Definition at line 97 of file forward_art_zernike3d.cpp.

98 {
99  addUsageLine("Template-based canonical volume reconstruction through Zernike3D coefficients");
100  defaultComments["-i"].clear();
101  defaultComments["-i"].addComment("Metadata with initial alignment");
102  defaultComments["-o"].clear();
103  defaultComments["-o"].addComment("Refined volume");
105  addParamsLine(" [--ref <volume=\"\">] : Reference volume");
106  addParamsLine(" [--mask <m=\"\">] : Mask reference volume");
107  addParamsLine(" [--odir <outputDir=\".\">] : Output directory");
108  addParamsLine(" [--sampling <Ts=1>] : Sampling rate (A/pixel)");
109  addParamsLine(" [--RDef <r=-1>] : Maximum radius of the deformation (px). -1=Half of volume size");
110  addParamsLine(" [--l1 <l1=3>] : Degree Zernike Polynomials=1,2,3,...");
111  addParamsLine(" [--l2 <l2=2>] : Harmonical depth of the deformation=1,2,3,...");
112  addParamsLine(" [--blobr <b=4>] : Blob radius for forward mapping splatting");
113  addParamsLine(" [--step <step=1>] : Voxel index step");
114  addParamsLine(" [--sigma <step=0.25>] : Gaussian sigma");
115  addParamsLine(" [--useZernike] : Correct heterogeneity with Zernike3D coefficients");
116  addParamsLine(" [--useCTF] : Correct CTF during ART reconstruction");
117  addParamsLine(" [--phaseFlipped] : Input images have been phase flipped");
118  addParamsLine(" [--regularization <l=0.01>] : ART regularization weight");
119  addParamsLine(" [--niter <n=1>] : Number of ART iterations");
120  addParamsLine(" [--save_iter <s=0>] : Save intermidiate volume after #save_iter iterations");
121  addParamsLine(" [--sort_last <N=2>] : The algorithm sorts projections in the most orthogonally possible way. ");
122  addParamsLine(" : The most orthogonal way is defined as choosing the projection which maximizes the ");
123  addParamsLine(" : dot product with the N previous inserted projections. Use -1 to sort with all ");
124  addParamsLine(" : previous projections");
125  addParamsLine(" [--resume] : Resume processing");
126  addExampleLine("A typical use is:", false);
127  addExampleLine("xmipp_forward_art_zernike3d -i anglesFromContinuousAssignment.xmd --ref reference.vol -o assigned_anglesAndDeformations.xmd --l1 3 --l2 2");
128 }
void addExampleLine(const char *example, bool verbatim=true)
void addUsageLine(const char *line, bool verbatim=false)
void addParamsLine(const String &line)
std::map< String, CommentList > defaultComments
Definition: xmipp_program.h:83

◆ deformVol()

void ProgForwardArtZernike3D::deformVol ( MultidimArray< double > &  mP,
MultidimArray< double > &  mW,
const MultidimArray< double > &  mV,
double  rot,
double  tilt,
double  psi 
)

Deform a volumen using Zernike-Spherical harmonic basis.

◆ fillVectorTerms()

void ProgForwardArtZernike3D::fillVectorTerms ( int  l1,
int  l2,
Matrix1D< int > &  vL1,
Matrix1D< int > &  vN,
Matrix1D< int > &  vL2,
Matrix1D< int > &  vM 
)

Zernike and SPH coefficients allocation.

Definition at line 331 of file forward_art_zernike3d.cpp.

333 {
334  int idx = 0;
335  vL1.initZeros(vecSize);
336  vN.initZeros(vecSize);
337  vL2.initZeros(vecSize);
338  vM.initZeros(vecSize);
339  for (int h = 0; h <= l2; h++)
340  {
341  int totalSPH = 2 * h + 1;
342  int aux = std::floor(totalSPH / 2);
343  for (int l = h; l <= l1; l += 2)
344  {
345  for (int m = 0; m < totalSPH; m++)
346  {
347  VEC_ELEM(vL1, idx) = l;
348  VEC_ELEM(vN, idx) = h;
349  VEC_ELEM(vL2, idx) = h;
350  VEC_ELEM(vM, idx) = m - aux;
351  idx++;
352  }
353  }
354  }
355 }
#define VEC_ELEM(v, i)
Definition: matrix1d.h:245
__host__ __device__ float2 floor(const float2 v)
void initZeros()
Definition: matrix1d.h:592
int m

◆ finishProcessing()

void ProgForwardArtZernike3D::finishProcessing ( )
virtual

Reimplemented from XmippMetadataProgram.

Definition at line 255 of file forward_art_zernike3d.cpp.

256 {
257  recoverVol();
258  Vout.write(fnVolO);
259 }
void write(const FileName &name="", size_t select_img=ALL_IMAGES, bool isStack=false, int mode=WRITE_OVERWRITE, CastWriteMode castMode=CW_CAST, int _swapWrite=0)

◆ numCoefficients()

void ProgForwardArtZernike3D::numCoefficients ( int  l1,
int  l2,
int &  vecSize 
)

Length of coefficients vector.

Definition at line 313 of file forward_art_zernike3d.cpp.

314 {
315  for (int h = 0; h <= l2; h++)
316  {
317  int numSPH = 2 * h + 1;
318  int count = l1 - h + 1;
319  int numEven = (count >> 1) + (count & 1 && !(h & 1));
320  if (h % 2 == 0)
321  {
322  vecSize += numSPH * numEven;
323  }
324  else
325  {
326  vecSize += numSPH * (l1 - h + 1 - numEven);
327  }
328  }
329 }

◆ preProcess()

void ProgForwardArtZernike3D::preProcess ( )
virtual

Produce side info. An exception is thrown if any of the files is not found

Reimplemented from XmippMetadataProgram.

Definition at line 130 of file forward_art_zernike3d.cpp.

131 {
132 
133  // Check that metadata has all information neede
134  if (!getInputMd()->containsLabel(MDL_ANGLE_ROT) ||
135  !getInputMd()->containsLabel(MDL_ANGLE_TILT) ||
136  !getInputMd()->containsLabel(MDL_ANGLE_PSI))
137  {
138  REPORT_ERROR(ERR_MD_MISSINGLABEL, "Input metadata projection angles are missing. Exiting...");
139  }
140 
141  if (fnVolR != "")
142  {
143  V.read(fnVolR);
144  }
145  else
146  {
147  FileName fn_first_image;
148  Image<double> first_image;
149  getInputMd()->getRow(1)->getValue(MDL_IMAGE, fn_first_image);
150  first_image.read(fn_first_image);
151  size_t Xdim_first = XSIZE(first_image());
152  V().initZeros(Xdim_first, Xdim_first, Xdim_first);
153  }
154  V().setXmippOrigin();
155 
156  Xdim = XSIZE(V());
157 
158  Vout().initZeros(V());
159  Vout().setXmippOrigin();
160 
161  if (resume && fnVolO.exists())
162  {
164  }
165  else
166  {
167  Vrefined() = V();
168  }
169  Vrefined().setXmippOrigin();
170 
171  if (RmaxDef < 0)
172  RmaxDef = Xdim / 2;
173 
174  // Transformation matrix
175  A.initIdentity(3);
176 
177  // CTF Filter
180  FilterCTF.ctf.enable_CTFnoise = false;
181  FilterCTF.ctf.enable_CTF = true;
182 
183  // Area where Zernike3D basis is computed (and volume is updated)
184  // Read Reference mask if avalaible (otherwise sphere of radius RmaxDef is used)
185  Mask mask;
186  mask.type = BINARY_CIRCULAR_MASK;
187  mask.mode = INNER_MASK;
188  if (fnMaskR != "")
189  {
190  Image<double> aux;
191  aux.read(fnMaskR);
192  typeCast(aux(), Vmask);
194  double Rmax2 = RmaxDef * RmaxDef;
195  for (int k = STARTINGZ(Vmask); k <= FINISHINGZ(Vmask); k++)
196  {
197  for (int i = STARTINGY(Vmask); i <= FINISHINGY(Vmask); i++)
198  {
199  for (int j = STARTINGX(Vmask); j <= FINISHINGX(Vmask); j++)
200  {
201  double r2 = k * k + i * i + j * j;
202  if (r2 >= Rmax2)
203  A3D_ELEM(Vmask, k, i, j) = 0;
204  }
205  }
206  }
207  }
208  else
209  {
210  mask.R1 = RmaxDef;
211  mask.generate_mask(V());
212  Vmask = mask.get_binary_mask();
214  }
215 
216  // Area Zernike3D in 2D
217  mask.R1 = RmaxDef;
218  mask.generate_mask(XSIZE(V()), XSIZE(V()));
219  mask2D = mask.get_binary_mask();
221 
222  vecSize = 0;
224  fillVectorTerms(L1, L2, vL1, vN, vL2, vM);
225 
226  // createWorkFiles();
227  initX = STARTINGX(Vrefined());
228  endX = FINISHINGX(Vrefined());
229  initY = STARTINGY(Vrefined());
230  endY = FINISHINGY(Vrefined());
231  initZ = STARTINGZ(Vrefined());
232  endZ = FINISHINGZ(Vrefined());
233 
236  filter.w1=sigma;
239  filter2.w1=sigma;
240 
241  // Blob
242  blob.radius = blob_r; // Blob radius in voxels
243  blob.order = 2; // Order of the Bessel function
244  blob.alpha = 3.6; // Smoothness parameter
245 
246  sigma4 = 2 * sigma;
253 }
Matrix1D< double > gaussianProjectionTable
Rotation angle of an image (double,degrees)
double alpha
Smoothness parameter.
Definition: blobs.h:121
CTFDescription ctf
#define FINISHINGX(v)
#define CTFINV
#define REPORT_ERROR(nerr, ErrormMsg)
Definition: xmipp_error.h:211
Definition: mask.h:360
void sqrt(Image< double > &op)
Tilting angle of an image (double,degrees)
Special label to be used when gathering MDs in MpiMetadataPrograms.
#define FINISHINGZ(v)
MultidimArray< int > Vmask
#define STARTINGX(v)
#define i
ql0001_ & k(htemp+1),(cvec+1),(atemp+1),(bj+1),(bl+1),(bu+1),(x+1),(clamda+1), &iout, infoqp, &zero,(w+1), &lenw,(iw+1), &leniw, &glob_grd.epsmac
Matrix1D< double > gaussianProjectionTable2
bool enable_CTF
Enable CTF part.
Definition: ctf.h:275
#define STARTINGY(v)
#define A3D_ELEM(V, k, i, j)
#define FOR_ALL_ELEMENTS_IN_MATRIX1D(v)
Definition: matrix1d.h:72
#define REALGAUSSIANZ
#define CEIL(x)
Definition: xmipp_macros.h:225
void fillVectorTerms(int l1, int l2, Matrix1D< int > &vL1, Matrix1D< int > &vN, Matrix1D< int > &vL2, Matrix1D< int > &vM)
Zernike and SPH coefficients allocation.
MultidimArray< int > mask2D
double R1
Definition: mask.h:413
void resize(size_t Xdim, bool copy=true)
Definition: matrix1d.h:410
#define XSIZE(v)
#define RAISED_COSINE
Missing expected label.
Definition: xmipp_error.h:158
int type
Definition: mask.h:402
virtual std::unique_ptr< MDRow > getRow(size_t id)=0
bool exists() const
void numCoefficients(int l1, int l2, int &vecSize)
Length of coefficients vector.
#define j
#define FINISHINGY(v)
void generate_mask(bool apply_geo=false)
Definition: mask.cpp:1577
void typeCast(const Matrix1D< T1 > &v1, Matrix1D< T2 > &v2)
Definition: matrix1d.h:1227
#define BINARY_CIRCULAR_MASK
Definition: mask.h:365
int order
Derivation order and Bessel function order.
Definition: blobs.h:118
int read(const FileName &name, DataMode datamode=DATA, size_t select_img=ALL_IMAGES, bool mapData=false, int mode=WRITE_READONLY)
float r2
bool enable_CTFnoise
Enable CTFnoise part.
Definition: ctf.h:273
const MultidimArray< int > & get_binary_mask() const
Definition: mask.h:707
#define REALGAUSSIANZ2
#define STARTINGZ(v)
double radius
Spatial radius in Universal System units.
Definition: blobs.h:115
Name of an image (std::string)
#define LOWPASS
void initIdentity()
Definition: matrix2d.h:673
int mode
Definition: mask.h:407
constexpr int INNER_MASK
Definition: mask.h:47
double gaussian1D(double x, double sigma, double mu)

◆ processImage()

void ProgForwardArtZernike3D::processImage ( const FileName fnImg,
const FileName fnImgOut,
const MDRow rowIn,
MDRow rowOut 
)
virtual

Predict angles and shift. At the input the pose parameters must have an initial guess of the parameters. At the output they have the estimated pose.

Implements XmippMetadataProgram.

Definition at line 262 of file forward_art_zernike3d.cpp.

263 {
264  flagEnabled = 1;
265 
266  int img_enabled;
267  rowIn.getValue(MDL_ENABLED, img_enabled);
268  if (img_enabled == -1) return;
269 
270  rowIn.getValue(MDL_ANGLE_ROT, rot);
271  rowIn.getValue(MDL_ANGLE_TILT, tilt);
272  rowIn.getValue(MDL_ANGLE_PSI, psi);
273  rowIn.getValueOrDefault(MDL_SHIFT_X, shiftX, 0.0);
274  rowIn.getValueOrDefault(MDL_SHIFT_Y, shiftY, 0.0);
275  std::vector<double> vectortemp;
276  if (useZernike)
277  {
278  rowIn.getValue(MDL_SPH_COEFFICIENTS, vectortemp);
279  std::vector<double> vec(vectortemp.begin(), vectortemp.end());
280  clnm = vec;
281  }
282  rowIn.getValueOrDefault(MDL_FLIP, flip, false);
283 
285  {
286  hasCTF = true;
287  FilterCTF.ctf.readFromMdRow(rowIn, false);
288  FilterCTF.ctf.Tm = Ts;
290  }
291  else
292  hasCTF = false;
293  MAT_ELEM(A, 0, 2) = shiftX;
294  MAT_ELEM(A, 1, 2) = shiftY;
295  MAT_ELEM(A, 0, 0) = 1;
296  MAT_ELEM(A, 0, 1) = 0;
297  MAT_ELEM(A, 1, 0) = 0;
298  MAT_ELEM(A, 1, 1) = 1;
299 
300  if (verbose >= 2)
301  std::cout << "Processing " << fnImg << std::endl;
302 
303  I.read(fnImg);
304  I().setXmippOrigin();
305 
306  // Forward Model
307  artModel<Direction::Forward>();
308 
309  // ART update
310  artModel<Direction::Backward>();
311 }
std::vector< double > clnm
Rotation angle of an image (double,degrees)
Defocus U (Angstroms)
CTFDescription ctf
Tilting angle of an image (double,degrees)
Shift for the image in the X axis (double)
Special label to be used when gathering MDs in MpiMetadataPrograms.
Name for the CTF Model (std::string)
Is this image enabled? (int [-1 or 1])
#define MAT_ELEM(m, i, j)
Definition: matrix2d.h:116
T & getValue(MDLabel label)
Flip the image? (bool)
double Tm
Sampling rate (A/pixel)
Definition: ctf.h:240
void readFromMdRow(const MDRow &row, bool disable_if_not_K=true)
Definition: ctf.cpp:1172
int verbose
Verbosity level.
const T & getValueOrDefault(MDLabel label, const T &def) const
virtual bool containsLabel(MDLabel label) const =0
void produceSideInfo()
Produce Side information.
Definition: ctf.cpp:1392
Deformation coefficients.
int read(const FileName &name, DataMode datamode=DATA, size_t select_img=ALL_IMAGES, bool mapData=false, int mode=WRITE_READONLY)
Shift for the image in the Y axis (double)

◆ readParams()

void ProgForwardArtZernike3D::readParams ( )
virtual

Read argument from command line.

Reimplemented from XmippMetadataProgram.

Definition at line 46 of file forward_art_zernike3d.cpp.

47 {
49  fnVolR = getParam("--ref");
50  fnMaskR = getParam("--mask");
51  fnOutDir = getParam("--odir");
52  RmaxDef = getIntParam("--RDef");
53  phaseFlipped = checkParam("--phaseFlipped");
54  useCTF = checkParam("--useCTF");
55  Ts = getDoubleParam("--sampling");
56  L1 = getIntParam("--l1");
57  L2 = getIntParam("--l2");
58  blob_r = getDoubleParam("--blobr");
59  loop_step = getIntParam("--step");
60  sigma = getDoubleParam("--sigma");
61  useZernike = checkParam("--useZernike");
62  lambda = getDoubleParam("--regularization");
63  resume = checkParam("--resume");
64  niter = getIntParam("--niter");
65  save_iter = getIntParam("--save_iter");
66  sort_last_N = getIntParam("--sort_last");
67  FileName outPath = getParam("-o");
68  outPath = outPath.afterLastOf("/");
69  fnVolO = fnOutDir + "/" + outPath;
70 }
double getDoubleParam(const char *param, int arg=0)
FileName afterLastOf(const String &str) const
FileName fnOutDir
Output directory.
const char * getParam(const char *param, int arg=0)
bool checkParam(const char *param)
int getIntParam(const char *param, int arg=0)

◆ recoverVol()

void ProgForwardArtZernike3D::recoverVol ( )

Definition at line 406 of file forward_art_zernike3d.cpp.

407 {
408  // Find the part of the volume that must be updated
409  auto &mVout = Vout();
410  const auto &mV = Vrefined();
411  mVout.initZeros(mV);
412  mVout = mV;
413 }

◆ show()

void ProgForwardArtZernike3D::show ( )

Show.

Definition at line 73 of file forward_art_zernike3d.cpp.

74 {
75  if (!verbose)
76  return;
78  std::cout
79  << "Output directory: " << fnOutDir << std::endl
80  << "Reference volume: " << fnVolR << std::endl
81  << "Reference mask: " << fnMaskR << std::endl
82  << "Sampling: " << Ts << std::endl
83  << "Max. Radius Deform. " << RmaxDef << std::endl
84  << "Zernike Degree: " << L1 << std::endl
85  << "SH Degree: " << L2 << std::endl
86  << "Blob radius: " << blob_r << std::endl
87  << "Step: " << loop_step << std::endl
88  << "Correct CTF: " << useCTF << std::endl
89  << "Correct heretogeneity: " << useZernike << std::endl
90  << "Phase flipped: " << phaseFlipped << std::endl
91  << "Regularization: " << lambda << std::endl
92  << "Number of iterations: " << niter << std::endl
93  << "Save every # iterations: " << save_iter << std::endl;
94 }
FileName fnOutDir
Output directory.
int verbose
Verbosity level.
void show() const override

Member Data Documentation

◆ A

Matrix2D<double> ProgForwardArtZernike3D::A

Definition at line 98 of file forward_art_zernike3d.h.

◆ blob

struct blobtype ProgForwardArtZernike3D::blob

Definition at line 132 of file forward_art_zernike3d.h.

◆ blob_r

double ProgForwardArtZernike3D::blob_r

Definition at line 133 of file forward_art_zernike3d.h.

◆ clnm

std::vector<double> ProgForwardArtZernike3D::clnm

Definition at line 116 of file forward_art_zernike3d.h.

◆ ctf

CTFDescription ProgForwardArtZernike3D::ctf

Definition at line 110 of file forward_art_zernike3d.h.

◆ current_image

size_t ProgForwardArtZernike3D::current_image

Definition at line 125 of file forward_art_zernike3d.h.

◆ current_iter

int ProgForwardArtZernike3D::current_iter

Definition at line 127 of file forward_art_zernike3d.h.

◆ current_save_iter

int ProgForwardArtZernike3D::current_save_iter

Definition at line 122 of file forward_art_zernike3d.h.

◆ defocusAngle

double ProgForwardArtZernike3D::defocusAngle

Definition at line 108 of file forward_art_zernike3d.h.

◆ defocusU

double ProgForwardArtZernike3D::defocusU

Definition at line 108 of file forward_art_zernike3d.h.

◆ defocusV

double ProgForwardArtZernike3D::defocusV

Definition at line 108 of file forward_art_zernike3d.h.

◆ endX

int ProgForwardArtZernike3D::endX

Definition at line 129 of file forward_art_zernike3d.h.

◆ endY

int ProgForwardArtZernike3D::endY

Definition at line 129 of file forward_art_zernike3d.h.

◆ endZ

int ProgForwardArtZernike3D::endZ

Definition at line 129 of file forward_art_zernike3d.h.

◆ filter

FourierFilter ProgForwardArtZernike3D::filter

Definition at line 142 of file forward_art_zernike3d.h.

◆ filter2

FourierFilter ProgForwardArtZernike3D::filter2

Definition at line 142 of file forward_art_zernike3d.h.

◆ FilterCTF

FourierFilter ProgForwardArtZernike3D::FilterCTF

Definition at line 112 of file forward_art_zernike3d.h.

◆ flagEnabled

int ProgForwardArtZernike3D::flagEnabled

Definition at line 72 of file forward_art_zernike3d.h.

◆ flip

bool ProgForwardArtZernike3D::flip

Definition at line 104 of file forward_art_zernike3d.h.

◆ fnMaskR

FileName ProgForwardArtZernike3D::fnMaskR

Filename of the reference volume mask

Definition at line 46 of file forward_art_zernike3d.h.

◆ fnOutDir

FileName ProgForwardArtZernike3D::fnOutDir

Output directory.

Definition at line 50 of file forward_art_zernike3d.h.

◆ fnVolO

FileName ProgForwardArtZernike3D::fnVolO

Filename of the refined volume

Definition at line 48 of file forward_art_zernike3d.h.

◆ fnVolR

FileName ProgForwardArtZernike3D::fnVolR

Filename of the reference volume

Definition at line 44 of file forward_art_zernike3d.h.

◆ gaussianProjectionTable

Matrix1D<double> ProgForwardArtZernike3D::gaussianProjectionTable

Definition at line 136 of file forward_art_zernike3d.h.

◆ gaussianProjectionTable2

Matrix1D<double> ProgForwardArtZernike3D::gaussianProjectionTable2

Definition at line 139 of file forward_art_zernike3d.h.

◆ hasCTF

bool ProgForwardArtZernike3D::hasCTF

Definition at line 106 of file forward_art_zernike3d.h.

◆ I

Image<double> ProgForwardArtZernike3D::I

Definition at line 88 of file forward_art_zernike3d.h.

◆ Idiff

Image<double> ProgForwardArtZernike3D::Idiff

Definition at line 96 of file forward_art_zernike3d.h.

◆ Ifilteredp

Image<double> ProgForwardArtZernike3D::Ifilteredp

Definition at line 86 of file forward_art_zernike3d.h.

◆ ignoreCTF

bool ProgForwardArtZernike3D::ignoreCTF

Definition at line 62 of file forward_art_zernike3d.h.

◆ initX

int ProgForwardArtZernike3D::initX

Definition at line 129 of file forward_art_zernike3d.h.

◆ initY

int ProgForwardArtZernike3D::initY

Definition at line 129 of file forward_art_zernike3d.h.

◆ initZ

int ProgForwardArtZernike3D::initZ

Definition at line 129 of file forward_art_zernike3d.h.

◆ L1

int ProgForwardArtZernike3D::L1

Degrees of Zernike polynomials and spherical harmonics

Definition at line 52 of file forward_art_zernike3d.h.

◆ L2

int ProgForwardArtZernike3D::L2

Definition at line 52 of file forward_art_zernike3d.h.

◆ lambda

double ProgForwardArtZernike3D::lambda

Definition at line 64 of file forward_art_zernike3d.h.

◆ loop_step

int ProgForwardArtZernike3D::loop_step

Definition at line 131 of file forward_art_zernike3d.h.

◆ mask2D

MultidimArray<int> ProgForwardArtZernike3D::mask2D

Definition at line 82 of file forward_art_zernike3d.h.

◆ niter

int ProgForwardArtZernike3D::niter

Definition at line 78 of file forward_art_zernike3d.h.

◆ num_images

size_t ProgForwardArtZernike3D::num_images

Definition at line 124 of file forward_art_zernike3d.h.

◆ ordered_list

MultidimArray<size_t> ProgForwardArtZernike3D::ordered_list

Definition at line 120 of file forward_art_zernike3d.h.

◆ P

Image<double> ProgForwardArtZernike3D::P

Definition at line 92 of file forward_art_zernike3d.h.

◆ phaseFlipped

bool ProgForwardArtZernike3D::phaseFlipped

Definition at line 60 of file forward_art_zernike3d.h.

◆ psi

double ProgForwardArtZernike3D::psi

Definition at line 100 of file forward_art_zernike3d.h.

◆ resume

bool ProgForwardArtZernike3D::resume

Resume computations

Definition at line 76 of file forward_art_zernike3d.h.

◆ RmaxDef

int ProgForwardArtZernike3D::RmaxDef

Maximum radius

Definition at line 58 of file forward_art_zernike3d.h.

◆ rot

double ProgForwardArtZernike3D::rot

Definition at line 100 of file forward_art_zernike3d.h.

◆ save_iter

int ProgForwardArtZernike3D::save_iter

Definition at line 66 of file forward_art_zernike3d.h.

◆ shiftX

double ProgForwardArtZernike3D::shiftX

Definition at line 102 of file forward_art_zernike3d.h.

◆ shiftY

double ProgForwardArtZernike3D::shiftY

Definition at line 102 of file forward_art_zernike3d.h.

◆ showOptimization

bool ProgForwardArtZernike3D::showOptimization

Definition at line 118 of file forward_art_zernike3d.h.

◆ sigma

double ProgForwardArtZernike3D::sigma

Definition at line 134 of file forward_art_zernike3d.h.

◆ sigma4

double ProgForwardArtZernike3D::sigma4

Definition at line 134 of file forward_art_zernike3d.h.

◆ sort_last_N

int ProgForwardArtZernike3D::sort_last_N

Definition at line 80 of file forward_art_zernike3d.h.

◆ tilt

double ProgForwardArtZernike3D::tilt

Definition at line 100 of file forward_art_zernike3d.h.

◆ Ts

double ProgForwardArtZernike3D::Ts

Sampling rate

Definition at line 56 of file forward_art_zernike3d.h.

◆ useCTF

bool ProgForwardArtZernike3D::useCTF

Definition at line 68 of file forward_art_zernike3d.h.

◆ useZernike

bool ProgForwardArtZernike3D::useZernike

Definition at line 70 of file forward_art_zernike3d.h.

◆ V

Image<double> ProgForwardArtZernike3D::V

Definition at line 86 of file forward_art_zernike3d.h.

◆ vecSize

int ProgForwardArtZernike3D::vecSize

Definition at line 114 of file forward_art_zernike3d.h.

◆ vL1

Matrix1D<int> ProgForwardArtZernike3D::vL1

Zernike and SPH coefficients vectors

Definition at line 54 of file forward_art_zernike3d.h.

◆ vL2

Matrix1D<int> ProgForwardArtZernike3D::vL2

Definition at line 54 of file forward_art_zernike3d.h.

◆ vM

Matrix1D<int> ProgForwardArtZernike3D::vM

Definition at line 54 of file forward_art_zernike3d.h.

◆ Vmask

MultidimArray<int> ProgForwardArtZernike3D::Vmask

Definition at line 90 of file forward_art_zernike3d.h.

◆ vN

Matrix1D<int> ProgForwardArtZernike3D::vN

Definition at line 54 of file forward_art_zernike3d.h.

◆ Vout

Image<double> ProgForwardArtZernike3D::Vout

Definition at line 86 of file forward_art_zernike3d.h.

◆ Vrefined

Image<double> ProgForwardArtZernike3D::Vrefined

Definition at line 86 of file forward_art_zernike3d.h.

◆ W

Image<double> ProgForwardArtZernike3D::W

Definition at line 94 of file forward_art_zernike3d.h.

◆ Xdim

size_t ProgForwardArtZernike3D::Xdim

Definition at line 84 of file forward_art_zernike3d.h.


The documentation for this class was generated from the following files: