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

#include <angular_project_library.h>

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

Public Member Functions

 ProgAngularProjectLibrary ()
 
 ProgAngularProjectLibrary (const ProgAngularProjectLibrary &)=delete
 
 ProgAngularProjectLibrary (const ProgAngularProjectLibrary &&)=delete
 
 ~ProgAngularProjectLibrary ()
 
ProgAngularProjectLibraryoperator= (const ProgAngularProjectLibrary &)=delete
 
ProgAngularProjectLibraryoperator= (const ProgAngularProjectLibrary &&)=delete
 
void readParams ()
 
void defineParams ()
 
void show ()
 
void run ()
 
void createGroupSamplingFiles (void)
 
void get_sym_vectors (std::vector< Matrix1D< double > > &sym_points)
 
void project_angle_vector (int my_init, int my_end, bool verbose=true)
 
- 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)
 
virtual void initComments ()
 
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 show () 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

Sampling mysampling
 
double sampling
 
double perturb_projection_vector
 
double psi_sampling
 
double max_tilt_angle
 
double min_tilt_angle
 
FileName output_file
 
FileName output_file_root
 
FileName fn_groups
 
projectionType projType
 Type of projection algorithm: More...
 
double paddFactor
 The padding factor for Fourier projection. More...
 
double maxFrequency
 The maximum frequency for Fourier projection. More...
 
int BSplineDeg
 The type of interpolation (NEAR. More...
 
FileName FnexperimentalImages
 
bool angular_distance_bool
 
bool remove_points_far_away_from_experimental_data_bool
 
double angular_distance
 
bool compute_closer_sampling_point_bool
 
bool compute_neighbors_bool
 
int symmetry
 
FileName fn_sym
 symmetry file for sampling More...
 
FileName fn_sym_neigh
 symmetry file for heighbors computation More...
 
int sym_order
 
FileName input_volume
 
Image< double > inputVol
 
int Xdim
 
int Ydim
 
bool only_winner
 
RealShearsInfoVshears
 
FourierProjectorVfourier
 
- Public Attributes inherited from XmippProgram
bool doRun
 
bool runWithoutArgs
 
int verbose
 Verbosity level. More...
 
int debug
 

Additional Inherited Members

- Protected Member Functions inherited from XmippProgram
void defineCommons ()
 
- 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

Parameter class for the projection library program

Definition at line 47 of file angular_project_library.h.

Constructor & Destructor Documentation

◆ ProgAngularProjectLibrary() [1/3]

ProgAngularProjectLibrary::ProgAngularProjectLibrary ( )

fil vector with symmetry axis Empty constructor

sampling object 1 by default

Definition at line 31 of file angular_project_library.cpp.

32 {
35  Vshears=nullptr;
36  Vfourier=nullptr;
37 
38 }
void setSampling(double sampling)
Definition: sampling.cpp:121

◆ ProgAngularProjectLibrary() [2/3]

ProgAngularProjectLibrary::ProgAngularProjectLibrary ( const ProgAngularProjectLibrary )
delete

◆ ProgAngularProjectLibrary() [3/3]

ProgAngularProjectLibrary::ProgAngularProjectLibrary ( const ProgAngularProjectLibrary &&  )
delete

◆ ~ProgAngularProjectLibrary()

ProgAngularProjectLibrary::~ProgAngularProjectLibrary ( )

Destructor

Definition at line 40 of file angular_project_library.cpp.

41 {
42  delete Vshears;
43  delete Vfourier;
44 }

Member Function Documentation

◆ createGroupSamplingFiles()

void ProgAngularProjectLibrary::createGroupSamplingFiles ( void  )

Create separate sampling files for subsets of the -experimental_images docfile

Definition at line 406 of file angular_project_library.cpp.

407 {
408 
409  //#define DEBUGTIME
410 #ifdef DEBUGTIME
411  time_t start,end;
412  double time_dif;
413  time (&start);
414 #endif
415 
416  //load txt file
418 #ifdef DEBUGTIME
419 
420  time (&end);
421  time_dif = difftime (end,start);
422  start=end;
423  printf ("re-read entire sampling file after %.2lf seconds\n", time_dif );
424 #endif
425 
426  StringVector blockList;
428  FileName fn_temp, fn_exp;
429  FileName my_output_file_root;
430  MetaDataVec SFBlock;
431 
433  int igrp=1;
434  for (StringVector::iterator it= blockList.begin();
435  it!=blockList.end(); it++,igrp++)
436  {
437  my_output_file_root.compose(output_file_root + "_group",igrp,"");
438  std::cerr<<"Writing group sampling file "<< my_output_file_root<<std::endl;
439 
440  fn_temp.compose(*it,fn_exp);
441  SFBlock.read(fn_temp);
442  if (SFBlock.size() > 0)//Do we really need this check?
443  //I guess so since user may have supplied a particular
444  //defocus classification. ROB
445  {
446  mysampling.fillExpDataProjectionDirectionByLR(fn_temp);//SFBlock@fn_groups
448  {
449  //find sampling point closer to experimental point (only 0) and bool
450  //and save docfile with this information
451  mysampling.findClosestSamplingPoint(fn_temp,my_output_file_root);
452  }
453 
454  //save saveSamplingFile
456  {
458  mysampling.saveSamplingFile(my_output_file_root,false);
459  }
460  }
461  }
462 #ifdef DEBUGTIME
463  time (&end);
464  time_dif = difftime (end,start);
465  start=end;
466  printf ("Written all group sampling files after %.2lf seconds\n", time_dif );
467 #endif
468 
469 
470 }
void read(const FileName &inFile, const std::vector< MDLabel > *desiredLabels=nullptr, bool decomposeStack=true) override
void findClosestSamplingPoint(const MetaData &DFi, const FileName &output_file_root)
Definition: sampling.cpp:1991
void computeNeighbors(bool only_winner=false)
Definition: sampling.cpp:1715
void compose(const String &str, const size_t no, const String &ext="")
void saveSamplingFile(const FileName &fn_base, bool write_vectors=true, bool write_sampling_sphere=false)
Definition: sampling.cpp:1495
void getBlocksInMetaDataFile(const FileName &inFile, StringVector &blockList)
std::vector< String > StringVector
Definition: xmipp_strings.h:35
size_t size() const override
void readSamplingFile(const FileName &infilename, bool read_vectors=true, bool read_sampling_sphere=false)
Definition: sampling.cpp:1592
FileName removeBlockName() const
void fillExpDataProjectionDirectionByLR(const MetaData &DFi)
Definition: sampling.cpp:2256

◆ defineParams()

void ProgAngularProjectLibrary::defineParams ( )
virtual

Usage message. This function shows the way of introducing these parameters.

Reimplemented from XmippProgram.

Definition at line 100 of file angular_project_library.cpp.

101 {
102  addUsageLine("Create a gallery of projections from a volume");
103  addUsageLine("+see [[http://xmipp.cnb.csic.es/twiki/bin/view/Xmipp/Angular_project_library_v3][here]] for further information about this program.");
104 
105  addParamsLine(" -i <input_volume_file> : Input Volume");
106  addParamsLine(" -o <output_file_name> : stack with output files");
107  addParamsLine(" [--sym <symmetry=c1>] : Symmetry to define sampling ");
108  addParamsLine(" : One of the 17 possible symmetries in");
109  addParamsLine(" : single particle electron microscopy");
110  addParamsLine(" [--sampling_rate <Ts=5>] : Distance in degrees between sampling points");
111  addParamsLine("==+Extra parameters==");
112  addParamsLine(" [--sym_neigh <symmetry>] : symmetry used to define neighbors, by default");
113  addParamsLine(" : same as sym");
114  addParamsLine(" [--psi_sampling <psi=360>] : sampling in psi, 360 -> no sampling in psi");
115  addParamsLine(" [--max_tilt_angle <tmax=180>] : maximum tilt angle in degrees");
116  addParamsLine(" [--min_tilt_angle <tmin=0>] : minimum tilt angle in degrees");
117  addParamsLine(" [--experimental_images <docfile=\"\">] : doc file with experimental data");
118  addParamsLine(" [--angular_distance <ang=20>] : Do not search a distance larger than...");
119  addParamsLine(" requires --experimental_images;");
120  addParamsLine(" [--closer_sampling_points] : create doc file with closest sampling points");
121  addParamsLine(" requires --experimental_images;");
122  addParamsLine(" [--near_exp_data] : remove points far away from experimental data");
123  addParamsLine(" requires --experimental_images;");
124  addParamsLine(" [--compute_neighbors] : create doc file with sampling point neighbors");
125  addParamsLine(" requires --angular_distance;");
126  addParamsLine(" [--method <method=fourier>] : Projection method");
127  addParamsLine(" where <method>");
128  addParamsLine(" real_space : Makes projections by ray tracing in real space");
129  addParamsLine(" fourier <pad=1> <maxfreq=0.25> <interp=bspline> : Takes a central slice in Fourier space");
130  addParamsLine(" : pad controls the padding factor, by default, the padded volume is");
131  addParamsLine(" : the same than the original volume. ");
132  addParamsLine(" : maxfreq is the maximum frequency for the pixels and by default ");
133  addParamsLine(" : pixels with frequency more than 0.25 are not considered.");
134  addParamsLine(" : interp is the method for interpolation and the values can be: ");
135  addParamsLine(" : nearest: Nearest Neighborhood ");
136  addParamsLine(" : linear: Linear ");
137  addParamsLine(" : bspline: Cubic BSpline ");
138  addParamsLine(" [--perturb <sigma=0.0>] : gaussian noise projection unit vectors ");
139  addParamsLine(" : a value=sin(sampling_rate)/4 ");
140  addParamsLine(" : may be a good starting point ");
141  addParamsLine(" [--groups <selfile=\"\">] : selfile with groups");
142  addParamsLine(" [--only_winner] : if set each experimental");
143  addParamsLine(" : point will have a unique neighbor");
144 
145  addExampleLine("Sample at 2 degrees and use c6 symmetry:", false);
146  addExampleLine("xmipp_angular_project_library -i in.vol -o out.stk --sym c6 --sampling_rate 2");
147 
148 }
void addExampleLine(const char *example, bool verbatim=true)
void addUsageLine(const char *line, bool verbatim=false)
void addParamsLine(const String &line)

◆ get_sym_vectors()

void ProgAngularProjectLibrary::get_sym_vectors ( std::vector< Matrix1D< double > > &  sym_points)

get all directions related by symmetry to (1,0,0)

◆ operator=() [1/2]

ProgAngularProjectLibrary& ProgAngularProjectLibrary::operator= ( const ProgAngularProjectLibrary )
delete

◆ operator=() [2/2]

ProgAngularProjectLibrary& ProgAngularProjectLibrary::operator= ( const ProgAngularProjectLibrary &&  )
delete

◆ project_angle_vector()

void ProgAngularProjectLibrary::project_angle_vector ( int  my_init,
int  my_end,
bool  verbose = true 
)

Project in all the directions between indexes init and end

Definition at line 194 of file angular_project_library.cpp.

195 {
196  Projection P;
197  FileName fn_proj;
198  int mySize;
199  int numberStepsPsi = 1;
200 
201  mySize=my_end-my_init+1;
202  if (psi_sampling < 360)
203  {
204  numberStepsPsi = (int) (359.99999/psi_sampling);
205  mySize *= numberStepsPsi;
206  }
207 
208  if (verbose)
209  init_progress_bar(mySize);
210  int myCounter=0;
211 
212  for (double mypsi=0;mypsi<360;mypsi += psi_sampling)
213  for (int i=0;i<my_init;i++)
214  myCounter++;
215 // if (shears && XSIZE(inputVol())!=0 && VShears==NULL)
216 // VShears=new RealShearsInfo(inputVol());
217  if (projType == SHEARS && XSIZE(inputVol())!=0 && Vshears==nullptr)
219  if (projType == FOURIER && XSIZE(inputVol())!=0 && Vfourier==nullptr)
221  paddFactor,
222  maxFrequency,
223  BSplineDeg);
224 
225  for (double mypsi=0;mypsi<360;mypsi += psi_sampling)
226  {
227  for (int i=my_init;i<=my_end;i++)
228  {
229  if (verbose)
230  progress_bar(i-my_init);
231 
235 
236 // if (shears)
237 // projectVolume(*VShears, P, Ydim, Xdim, rot,tilt,psi);
238 // else
239 // projectVolume(inputVol(), P, Ydim, Xdim, rot,tilt,psi);
240  if (projType == SHEARS)
241  projectVolume(*Vshears, P, Ydim, Xdim, rot, tilt, psi);
242  else if (projType == FOURIER)
243  projectVolume(*Vfourier, P, Ydim, Xdim, rot, tilt, psi);
244  else if (projType == REALSPACE)
245  projectVolume(inputVol(), P, Ydim, Xdim, rot, tilt, psi);
246 
247 
248  P.setEulerAngles(rot,tilt,psi);
250  P.write(output_file,(size_t) (numberStepsPsi * i + mypsi +1),true,WRITE_REPLACE);
251  }
252  }
253  if (verbose)
254  progress_bar(mySize);
255 }
void init_progress_bar(long total)
constexpr int FOURIER
int BSplineDeg
The type of interpolation (NEAR.
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)
Definition: project.h:42
#define i
#define XX(v)
Definition: matrix1d.h:85
double paddFactor
The padding factor for Fourier projection.
#define XSIZE(v)
void progress_bar(long rlen)
Structure for holding a volume.
void projectVolume(FourierProjector &projector, Projection &P, int Ydim, int Xdim, double rot, double tilt, double psi, const MultidimArray< double > *ctf)
int verbose
Verbosity level.
void setDataMode(DataMode mode)
#define YY(v)
Definition: matrix1d.h:93
void setEulerAngles(double rot, double tilt, double psi, const size_t n=0)
std::vector< Matrix1D< double > > no_redundant_sampling_points_angles
Definition: sampling.h:121
double psi(const double x)
projectionType projType
Type of projection algorithm:
#define ZZ(v)
Definition: matrix1d.h:101
double maxFrequency
The maximum frequency for Fourier projection.

◆ readParams()

void ProgAngularProjectLibrary::readParams ( )
virtual

Read from a command line. An exception might be thrown by any of the internal conversions, this would mean that there is an error in the command line and you might show a usage message.

Reimplemented from XmippProgram.

Definition at line 47 of file angular_project_library.cpp.

48 {
49  input_volume = getParam("-i");
50  output_file = getParam("-o");
52  fn_sym = getParam("--sym");
53  fn_sym_neigh=checkParam("--sym_neigh")?getParam("--sym_neigh"):fn_sym;
54  sampling = getDoubleParam("--sampling_rate");
55  psi_sampling = getDoubleParam("--psi_sampling");
56  max_tilt_angle = getDoubleParam("--max_tilt_angle");
57  min_tilt_angle = getDoubleParam("--min_tilt_angle");
58  angular_distance_bool = checkParam("--angular_distance");
61  {
62  FnexperimentalImages = getParam("--experimental_images");
63  angular_distance = getDoubleParam("--angular_distance");
64  }
65  compute_closer_sampling_point_bool= checkParam("--closer_sampling_points");
67  FnexperimentalImages = getParam("--experimental_images");
68  if (STR_EQUAL(getParam("--method"), "real_space"))
70  if (STR_EQUAL(getParam("--method"), "shears"))
71  projType = SHEARS;
72  if (STR_EQUAL(getParam("--method"), "fourier"))
73  {
74  projType = FOURIER;
75  paddFactor = getDoubleParam("--method", 1);
76  maxFrequency = getDoubleParam("--method", 2);
77  String degree = getParam("--method", 3);
78  if (degree == "nearest")
80  else if (degree == "linear")
82  else if (degree == "bspline")
84  else
85  REPORT_ERROR(ERR_ARG_BADCMDLINE, "The interpolation kernel can be : nearest, linear, bspline");
86  }
87 
88  //NOTE perturb in computed after the even sampling is computes
89  // and max tilt min tilt applied
91  compute_neighbors_bool=checkParam("--compute_neighbors");
94  FnexperimentalImages = getParam("--experimental_images");
95  fn_groups = getParam("--groups");
96  only_winner = checkParam("--only_winner");
97 }
Errors on command line parameters.
Definition: xmipp_error.h:112
constexpr int FOURIER
double getDoubleParam(const char *param, int arg=0)
#define REPORT_ERROR(nerr, ErrormMsg)
Definition: xmipp_error.h:211
int BSplineDeg
The type of interpolation (NEAR.
FileName fn_sym
symmetry file for sampling
Definition: project.h:42
FileName fn_sym_neigh
symmetry file for heighbors computation
const char * getParam(const char *param, int arg=0)
double paddFactor
The padding factor for Fourier projection.
#define STR_EQUAL(str1, str2)
Definition: xmipp_strings.h:42
FileName withoutExtension() const
std::string String
Definition: xmipp_strings.h:34
projectionType projType
Type of projection algorithm:
bool checkParam(const char *param)
double maxFrequency
The maximum frequency for Fourier projection.

◆ run()

void ProgAngularProjectLibrary::run ( )
virtual

Run.

Reimplemented from XmippProgram.

Definition at line 258 of file angular_project_library.cpp.

259 {
261  // PreRun for all nodes but not for all works
263  //only rank 0
265  show();
266  //all ranks
268  srand ( time(nullptr) );
269  //process the symmetry file
270  //only checks symmetry and set pg_order and pg_group, no memory allocation
273  (std::string)"Invalid symmetry" + fn_sym);
275  {
276  int my_seed;
277  my_seed=rand();
278  // set noise deviation and seed
280  }
281  if(angular_distance_bool!=0)
283  //true -> half_sphere
285  //only rank 0
286  //mysampling.createSymFile(fn_sym,symmetry, sym_order);
287  //all nodes
289  //store symmetry matrices, this is faster than computing them each time
291  //mpi_barrier here
292  //all working nodes must read symmetry file
293  //and experimental docfile if apropiate
294  //symmetry_file = symmetry + ".sym";
295  //SL.readSymmetryFile(symmetry_file)
296  // We first sample The whole sphere
297  // Then we remove point redundant due to sampling symmetry
298  // use old symmetry, this is geometric does not use L_R
300 
301  //=========================
302  //======================
303  //recompute symmetry with neigh symmetry
304  // If uncomment neighbour are not OK. BE CAREFUL
305 #define BREAKSIMMETRY
306 #ifdef BREAKSIMMETRY
309  (std::string)"Invalid neig symmetry" + fn_sym_neigh);
312 #endif
313 #undef BREAKSIMMETRY
314  //precompute product between symmetry matrices and experimental data
315  if (FnexperimentalImages.size() > 0)
317 
318  //remove points not close to experimental points, only for no symmetric cases
319  if (FnexperimentalImages.size() > 0 &&
321  {
322  // here we remove points no close to experimental data, neight symmetry must be use
324  }
326  {
327  //find sampling point closer to experimental point (only 0) and bool
328  //and save docfile with this information
329  // use neight symmetry
331  }
332  //only rank 0
333  //write docfile with vectors and angles
335  //all nodes
336  //If there is no reference available exit
338  inputVol().setXmippOrigin();
339  Xdim = XSIZE(inputVol());
340  Ydim = YSIZE(inputVol());
341 
343  {
344  // new symmetry
347  }
348  //release some memory
350 
351  unlink(output_file.c_str());
352 
353  //mpi master should divide doc in chuncks
354  //in this serial program there is a unique chunck
355  //angle information is in
356  //mysampling.no_redundant_sampling_points_vector[i]
357  //Run for all works
360 
361  //only rank 0 create sel file
362  MetaDataVec mySFin, mySFout;
363  FileName fn_temp;
364  mySFin.read(output_file_root+"_angles.doc");
365  size_t myCounter=0;
366  size_t id;
367  int ref;
368 
369  for (double mypsi=0;mypsi<360;mypsi += psi_sampling)
370  {
371  for (size_t objId : mySFin.ids())
372  {
373  double x,y,z, rot, tilt, psi;
374  mySFin.getValue(MDL_ANGLE_ROT,rot,objId);
375  mySFin.getValue(MDL_ANGLE_TILT,tilt,objId);
376  mySFin.getValue(MDL_ANGLE_PSI,psi,objId);
377  mySFin.getValue(MDL_X,x,objId);
378  mySFin.getValue(MDL_Y,y,objId);
379  mySFin.getValue(MDL_Z,z,objId);
380  mySFin.getValue(MDL_REF,ref,objId);
381  fn_temp.compose( ++myCounter,output_file);
382  id = mySFout.addObject();
383  mySFout.setValue(MDL_IMAGE,fn_temp,id);
384  mySFout.setValue(MDL_ENABLED,1,id);
385  mySFout.setValue(MDL_ANGLE_ROT,rot,id);
386  mySFout.setValue(MDL_ANGLE_TILT,tilt,id);
387  mySFout.setValue(MDL_ANGLE_PSI,psi+mypsi,id);
388  mySFout.setValue(MDL_X,x,id);
389  mySFout.setValue(MDL_Y,y,id);
390  mySFout.setValue(MDL_Z,z,id);
391  mySFout.setValue(MDL_SCALE,1.0,id);
392  mySFout.setValue(MDL_REF,ref,id);
393  }
394  }
395  mySFout.setComment("x,y,z refer to the coordinates of the unitary vector at direction given by the euler angles");
396  if (mySFout.size()>0)
397  mySFout.write(output_file_root+".doc");
398  else
399  std::cout << "There are no projections within the specified angular range and sampling" << std::endl;
400  unlink((output_file_root+"_angles.doc").c_str());
401 
402  if (fn_groups!="")
404 }
Rotation angle of an image (double,degrees)
#define YSIZE(v)
void setSampling(double sampling)
Definition: sampling.cpp:121
void removeRedundantPoints(const int symmetry, int sym_order)
Definition: sampling.cpp:691
void setNoise(double deviation, int my_seed=-1)
Definition: sampling.cpp:134
#define REPORT_ERROR(nerr, ErrormMsg)
Definition: xmipp_error.h:211
void createAsymUnitFile(const FileName &docfilename)
Definition: sampling.cpp:1440
Z component (double)
void setNeighborhoodRadius(double neighborhood)
Definition: sampling.cpp:140
Tilting angle of an image (double,degrees)
static double * y
void findClosestSamplingPoint(const MetaData &DFi, const FileName &output_file_root)
Definition: sampling.cpp:1991
int verbose
Definition: sampling.h:126
int readSymmetryFile(FileName fn_sym, double accuracy=SYM_ACCURACY)
Definition: symmetries.cpp:33
bool isSymmetryGroup(FileName fn_sym, int &pgGroup, int &pgOrder)
Definition: symmetries.cpp:601
FileName fn_sym
symmetry file for sampling
void computeNeighbors(bool only_winner=false)
Definition: sampling.cpp:1715
void compose(const String &str, const size_t no, const String &ext="")
Special label to be used when gathering MDs in MpiMetadataPrograms.
void saveSamplingFile(const FileName &fn_base, bool write_vectors=true, bool write_sampling_sphere=false)
Definition: sampling.cpp:1495
doublereal * x
Is this image enabled? (int [-1 or 1])
FileName fn_sym_neigh
symmetry file for heighbors computation
std::vector< Matrix1D< double > > exp_data_projection_direction_by_L_R
Definition: sampling.h:108
void computeSamplingPoints(bool only_half_sphere=true, double max_tilt=180, double min_tilt=0)
Definition: sampling.cpp:155
void project_angle_vector(int my_init, int my_end, bool verbose=true)
#define XSIZE(v)
double z
scaling factor for an image or volume (double)
int verbose
Verbosity level.
void fillLRRepository(void)
Definition: sampling.cpp:2216
Class to which the image belongs (int)
std::vector< Matrix1D< double > > no_redundant_sampling_points_angles
Definition: sampling.h:121
Y component (double)
double psi(const double x)
X component (double)
int read(const FileName &name, DataMode datamode=DATA, size_t select_img=ALL_IMAGES, bool mapData=false, int mode=WRITE_READONLY)
void fillExpDataProjectionDirectionByLR(const MetaData &DFi)
Definition: sampling.cpp:2256
void removePointsFarAwayFromExperimentalData()
Definition: sampling.cpp:1928
Incorrect value received.
Definition: xmipp_error.h:195
Name of an image (std::string)
SymList SL
Definition: sampling.h:138

◆ show()

void ProgAngularProjectLibrary::show ( )

Show parameters.

Definition at line 151 of file angular_project_library.cpp.

152 {
153  if (!verbose)
154  return;
155  std::cout << "output input_volume root: " << input_volume << std::endl
156  << "output files root: " << output_file_root << std::endl
157  << "Sampling rate: " << sampling << std::endl
158  << "symmetry sampling group: " << fn_sym << std::endl
159  << "symmetry neighbohound group: " << fn_sym_neigh << std::endl
160  << "max_tilt_angle: " << max_tilt_angle << std::endl
161  << "min_tilt_angle: " << min_tilt_angle << std::endl
162  << "psi_sampling: " << psi_sampling << std::endl
163  << "compute_neighbors: " << compute_neighbors_bool << std::endl
164  << "only_winner: " << only_winner << std::endl
165  << "verbose: " << verbose << std::endl
166  << "projection method: ";
167 
168  if (projType == FOURIER)
169  {
170  std::cout << " fourier " <<std::endl;
171  std::cout << " pad factor: " << paddFactor <<std::endl;
172  std::cout << " maxFrequency: " << maxFrequency <<std::endl;
173  std::cout << " interpolator: ";
175  std::cout << " nearest" <<std::endl;
177  std::cout << " linear" <<std::endl;
179  std::cout << " bspline" <<std::endl;
180  }
181  else if (projType == REALSPACE)
182  std::cout << " realspace " <<std::endl;
183 
185  std::cout << "angular_distance: " << angular_distance << std::endl;
186  if (FnexperimentalImages.size() > 0)
187  std::cout << "experimental_images: " << FnexperimentalImages << std::endl;
188  std::cout << "compute_closer_sampling_point_bool:" << compute_closer_sampling_point_bool
189  << std::endl;
191  std::cout << "perturb_projection_vector: " << perturb_projection_vector << std::endl;
192 }
constexpr int FOURIER
int BSplineDeg
The type of interpolation (NEAR.
FileName fn_sym
symmetry file for sampling
FileName fn_sym_neigh
symmetry file for heighbors computation
double paddFactor
The padding factor for Fourier projection.
int verbose
Verbosity level.
projectionType projType
Type of projection algorithm:
double maxFrequency
The maximum frequency for Fourier projection.

Member Data Documentation

◆ angular_distance

double ProgAngularProjectLibrary::angular_distance

enabled angular_distance

Definition at line 104 of file angular_project_library.h.

◆ angular_distance_bool

bool ProgAngularProjectLibrary::angular_distance_bool

enabled angular_distance

Definition at line 98 of file angular_project_library.h.

◆ BSplineDeg

int ProgAngularProjectLibrary::BSplineDeg

The type of interpolation (NEAR.

Definition at line 84 of file angular_project_library.h.

◆ compute_closer_sampling_point_bool

bool ProgAngularProjectLibrary::compute_closer_sampling_point_bool

enabled angular_distance

Definition at line 107 of file angular_project_library.h.

◆ compute_neighbors_bool

bool ProgAngularProjectLibrary::compute_neighbors_bool

enable is neighbors must be computed

Definition at line 110 of file angular_project_library.h.

◆ fn_groups

FileName ProgAngularProjectLibrary::fn_groups

Name of selfile for groups

Definition at line 75 of file angular_project_library.h.

◆ fn_sym

FileName ProgAngularProjectLibrary::fn_sym

symmetry file for sampling

Definition at line 118 of file angular_project_library.h.

◆ fn_sym_neigh

FileName ProgAngularProjectLibrary::fn_sym_neigh

symmetry file for heighbors computation

Definition at line 121 of file angular_project_library.h.

◆ FnexperimentalImages

FileName ProgAngularProjectLibrary::FnexperimentalImages

filename with experimental images angles. This information is used to generate only angles close to experimental images

Definition at line 95 of file angular_project_library.h.

◆ input_volume

FileName ProgAngularProjectLibrary::input_volume

Input volume name

Definition at line 127 of file angular_project_library.h.

◆ inputVol

Image<double> ProgAngularProjectLibrary::inputVol

volume to be projecte

Definition at line 130 of file angular_project_library.h.

◆ max_tilt_angle

double ProgAngularProjectLibrary::max_tilt_angle

maximun tilt angle

Definition at line 63 of file angular_project_library.h.

◆ maxFrequency

double ProgAngularProjectLibrary::maxFrequency

The maximum frequency for Fourier projection.

Definition at line 82 of file angular_project_library.h.

◆ min_tilt_angle

double ProgAngularProjectLibrary::min_tilt_angle

minimum tilt angle

Definition at line 66 of file angular_project_library.h.

◆ mysampling

Sampling ProgAngularProjectLibrary::mysampling

sampling object

Definition at line 51 of file angular_project_library.h.

◆ only_winner

bool ProgAngularProjectLibrary::only_winner

If true there will be only one neighbour per sampling point, the closest

Definition at line 139 of file angular_project_library.h.

◆ output_file

FileName ProgAngularProjectLibrary::output_file

output file

Definition at line 69 of file angular_project_library.h.

◆ output_file_root

FileName ProgAngularProjectLibrary::output_file_root

root for output files

Definition at line 72 of file angular_project_library.h.

◆ paddFactor

double ProgAngularProjectLibrary::paddFactor

The padding factor for Fourier projection.

Definition at line 80 of file angular_project_library.h.

◆ perturb_projection_vector

double ProgAngularProjectLibrary::perturb_projection_vector

Perturb angles of reference projections

Definition at line 57 of file angular_project_library.h.

◆ projType

projectionType ProgAngularProjectLibrary::projType

Type of projection algorithm:

Definition at line 78 of file angular_project_library.h.

◆ psi_sampling

double ProgAngularProjectLibrary::psi_sampling

produce projections in psi?

Definition at line 60 of file angular_project_library.h.

◆ remove_points_far_away_from_experimental_data_bool

bool ProgAngularProjectLibrary::remove_points_far_away_from_experimental_data_bool

remove points far away from experimental data

Definition at line 101 of file angular_project_library.h.

◆ sampling

double ProgAngularProjectLibrary::sampling

Sampling rate. Distance between sampling points in degrees

Definition at line 54 of file angular_project_library.h.

◆ sym_order

int ProgAngularProjectLibrary::sym_order

For infinite groups symmetry order

Definition at line 124 of file angular_project_library.h.

◆ symmetry

int ProgAngularProjectLibrary::symmetry

Symmetry. One of the 17 possible symmetries in single particle electron microscopy.

Definition at line 115 of file angular_project_library.h.

◆ Vfourier

FourierProjector* ProgAngularProjectLibrary::Vfourier

Definition at line 145 of file angular_project_library.h.

◆ Vshears

RealShearsInfo* ProgAngularProjectLibrary::Vshears

Definition at line 142 of file angular_project_library.h.

◆ Xdim

int ProgAngularProjectLibrary::Xdim

projection x dim

Definition at line 133 of file angular_project_library.h.

◆ Ydim

int ProgAngularProjectLibrary::Ydim

projection y dim

Definition at line 136 of file angular_project_library.h.


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