Xmipp  v3.23.11-Nereus
Public Member Functions | Public Attributes | List of all members
ProgVolumeValidationPCA Class Reference
Libraries » Reconstruction library » Validation of volume consistency with respect to the provided classes

#include <volume_validate_pca.h>

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

Public Member Functions

void readParams ()
 Read arguments from command line. More...
 
void defineParams ()
 Read arguments from command line. More...
 
void show ()
 
void run ()
 
void produceSideinfo ()
 Produce side info: fill arrays with relevant transformation matrices. More...
 
void reconstruct ()
 Reconstruct volume with all projections. More...
 
void modifyAngles ()
 Modified angles. More...
 
void evaluate ()
 Modified angles. 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)
 
ParamDef * getParamDef (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

FileName fnClasses
 
FileName fnSym
 
FileName fnAngles
 
int NVols
 
int NClasses
 
bool doProjMatch
 
MetaDataVec mdClasses
 
MetaDataVec mdReconstruction
 
- 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
 
ProgramDef * progDef
 Program definition and arguments parser. More...
 
std::map< String, CommentList > defaultComments
 
int argc
 Original command line arguments. More...
 
const char ** argv
 

Detailed Description

PCA validation volume parameters.

Definition at line 43 of file volume_validate_pca.h.

Member Function Documentation

◆ defineParams()

void ProgVolumeValidationPCA::defineParams ( )
virtual

Read arguments from command line.

Reimplemented from XmippProgram.

Definition at line 32 of file volume_validate_pca.cpp.

33 {
34  //usage
35  addUsageLine("Validate an obtained volume from a set of class averages");
36  //params
37  addParamsLine(" -i <md_file> : Metadata file with input classes");
38  addParamsLine(" [ --vol <file=\"\">] : Input volume");
39  addParamsLine(" [--sym <symfile=c1>] : Enforce symmetry in projections");
40  addParamsLine(" [--numVols <N=5>] : Number of intermediate volumes to generate");
41  addParamsLine(" [--numClasses <N=8>] : Number of classes to generate the intermediate volumes");
42 }
XmippProgram::addUsageLine
void addUsageLine(const char *line, bool verbatim=false)
Definition: xmipp_program.cpp:282
XmippProgram::addParamsLine
void addParamsLine(const String &line)
Definition: xmipp_program.cpp:305

◆ evaluate()

void ProgVolumeValidationPCA::evaluate ( )

Modified angles.

Definition at line 168 of file volume_validate_pca.cpp.

169 {
170  MetaDataVec md, md2;
171  FileName imag;
172  double rot, tilt, psi, rot2, tilt2, psi2;
173  std::vector<double> error;
174  size_t id, id2;
175 
176  md.read(fnClasses);
177  std::stringstream ss;
178 
179  MDMultiQuery query;
180  MetaDataVec auxMetadata;
181  double distance = 0;
182 
183  for (const auto& row : md)
184  {
185  row.getValue(MDL_IMAGE, imag);
186  row.getValue(MDL_ANGLE_ROT, rot);
187  row.getValue(MDL_ANGLE_TILT,tilt);
188  row.getValue(MDL_ANGLE_PSI,psi);
189 
190  for (int index=0; index<NVols;index++)
191  {
192  ss << index;
193  fnAngles=fnClasses.removeAllExtensions()+ss.str();
194  fnAngles+=".xmd";
195  ss.str(std::string());
196  md2.read(fnAngles);
197  MDValueEQ eq(MDL_IMAGE, imag);
198  query.addAndQuery(eq);
199  auxMetadata.importObjects(md2, eq);
200 
201  for (const auto& rowAux : auxMetadata)
202  {
203  rowAux.getValue(MDL_ANGLE_ROT, rot2);
204  rowAux.getValue(MDL_ANGLE_TILT,tilt2);
205  rowAux.getValue(MDL_ANGLE_PSI,psi2);
206  distance = Euler_distanceBetweenAngleSets(rot,tilt, psi,rot2, tilt2, psi2, true);
207  error.push_back(distance);
208  }
209  }
210  }
211 
212  double sum = std::accumulate(error.begin(), error.end(), 0.0);
213  double mean = sum / error.size();
214  double sq_sum = std::inner_product(error.begin(), error.end(), error.begin(), 0.0);
215  double stdev = std::sqrt(sq_sum / error.size() - mean * mean);
216 
217  std::cout << std::endl;
218  std::cout << "stdev of errors : " << stdev << std::endl;
219  std::cout << "mean of errors : " << mean << std::endl;
220 
221 }
MDL_ANGLE_ROT
Rotation angle of an image (double,degrees)
Definition: metadata_label.h:55
MDMultiQuery::addAndQuery
void addAndQuery(MDQuery &query)
Definition: metadata_query.h:344
MetaDataVec::read
void read(const FileName &inFile, const std::vector< MDLabel > *desiredLabels=nullptr, bool decomposeStack=true) override
Definition: metadata_vec.cpp:120
sqrt
void sqrt(Image< double > &op)
Definition: image_operate.cpp:210
MDL_ANGLE_TILT
Tilting angle of an image (double,degrees)
Definition: metadata_label.h:59
MDL_ANGLE_PSI
Special label to be used when gathering MDs in MpiMetadataPrograms.
Definition: metadata_label.h:51
FileName::removeAllExtensions
FileName removeAllExtensions() const
Definition: xmipp_filename.cpp:382
index
viol index
Definition: numerical_recipes.cpp:5197
MetaDataVec::importObjects
void importObjects(const MetaData &md, const std::vector< size_t > &objectsToAdd, bool doClear=true) override
Definition: metadata_vec.cpp:458
eq
void eq(Image< double > &op1, const Image< double > &op2)
Be careful with integer images for relational operations...due to double comparisons.
Definition: image_operate.cpp:149
Euler_distanceBetweenAngleSets
T Euler_distanceBetweenAngleSets(T rot1, T tilt1, T psi1, T rot2, T tilt2, T psi2, bool only_projdir)
Definition: geometry.cpp:681
distance
TYPE distance(struct Point_T *p, struct Point_T *q)
Definition: point.cpp:28
error
void error(char *s)
Definition: tools.cpp:107
alglib::psi
double psi(const double x)
Definition: specialfunctions.cpp:2492
MDL_IMAGE
Name of an image (std::string)
Definition: metadata_label.h:246

◆ modifyAngles()

void ProgVolumeValidationPCA::modifyAngles ( )

Modified angles.

Definition at line 74 of file volume_validate_pca.cpp.

75 {
76  //String args=formatString("-i %s -o %s --operate random_subset %d --mode overwrite",fnClasses.c_str(),fnAngles.c_str(),NClasses);
77  String args=formatString("-i %s -o %s --operate bootstrap --mode overwrite",fnClasses.c_str(),fnAngles.c_str());
78  String cmd=(String)"xmipp_metadata_utilities "+args;
79  if (system(cmd.c_str())==-1)
80  REPORT_ERROR(ERR_UNCLASSIFIED,"Cannot open shell");
81 
82  //Metadata with the well sampled projection and random projections assigned
83  MetaDataVec md;
84  double rot, tilt, psi;
85  size_t id;
86 
87  md.read(fnAngles);
88  randomize_random_generator();
89 
90  int anglePertur = 30;
91  for (auto& row : md)
92  {
93  row.getValue(MDL_ANGLE_ROT, rot);
94  row.getValue(MDL_ANGLE_TILT,tilt);
95  row.getValue(MDL_ANGLE_PSI,psi);
96 
97  row.setValue(MDL_ANGLE_ROT, (rnd_unif(-anglePertur,anglePertur))/2+rot);
98  row.setValue(MDL_ANGLE_TILT,(rnd_unif(-anglePertur,anglePertur))/2+tilt);
99  row.setValue(MDL_ANGLE_PSI, (rnd_unif(-anglePertur,anglePertur))/2+psi);
100 
101  md.setRow(row, row.id());
102  }
103 
104  md.write(fnAngles);
105 }
ERR_UNCLASSIFIED
Just to locate unclassified errors.
Definition: xmipp_error.h:192
MDL_ANGLE_ROT
Rotation angle of an image (double,degrees)
Definition: metadata_label.h:55
MetaDataVec::read
void read(const FileName &inFile, const std::vector< MDLabel > *desiredLabels=nullptr, bool decomposeStack=true) override
Definition: metadata_vec.cpp:120
REPORT_ERROR
#define REPORT_ERROR(nerr, ErrormMsg)
Definition: xmipp_error.h:211
MDL_ANGLE_TILT
Tilting angle of an image (double,degrees)
Definition: metadata_label.h:59
MDL_ANGLE_PSI
Special label to be used when gathering MDs in MpiMetadataPrograms.
Definition: metadata_label.h:51
rnd_unif
double rnd_unif()
Definition: xmipp_funcs.cpp:440
String
std::string String
Definition: xmipp_strings.h:34
alglib::psi
double psi(const double x)
Definition: specialfunctions.cpp:2492
formatString
String formatString(const char *format,...)
Definition: xmipp_strings.cpp:602
randomize_random_generator
unsigned int randomize_random_generator()
Definition: xmipp_funcs.cpp:407

◆ produceSideinfo()

void ProgVolumeValidationPCA::produceSideinfo ( )

Produce side info: fill arrays with relevant transformation matrices.

Definition at line 66 of file volume_validate_pca.cpp.

67 {
68  doProjMatch = false;
69  mdClasses.read(fnClasses);
70  mdClasses.removeDisabled();
71  getImageSize(mdClasses,xdim, ydim, zdim, ndim);
72 }
MetaDataVec::read
void read(const FileName &inFile, const std::vector< MDLabel > *desiredLabels=nullptr, bool decomposeStack=true) override
Definition: metadata_vec.cpp:120
getImageSize
void getImageSize(const MetaData &md, size_t &Xdim, size_t &Ydim, size_t &Zdim, size_t &Ndim, MDLabel image_label)
Definition: metadata_extension.cpp:273
MetaData::removeDisabled
virtual void removeDisabled()
Definition: metadata_base.cpp:123

◆ readParams()

void ProgVolumeValidationPCA::readParams ( )
virtual

Read arguments from command line.

Reimplemented from XmippProgram.

Definition at line 45 of file volume_validate_pca.cpp.

46 {
47  fnClasses = getParam("-i");
48  fnSym = getParam("--sym");
49  NVols = getIntParam("--numVols");
50  NClasses = getIntParam("--numClasses");
51 }
XmippProgram::getParam
const char * getParam(const char *param, int arg=0)
Definition: xmipp_program.cpp:321
XmippProgram::getIntParam
int getIntParam(const char *param, int arg=0)
Definition: xmipp_program.cpp:331

◆ reconstruct()

void ProgVolumeValidationPCA::reconstruct ( )

Reconstruct volume with all projections.

Definition at line 107 of file volume_validate_pca.cpp.

108 {
109  //number of threads in the reconstruction
110  int Nthr = 4;
111  float AngularSampling = 30;
112  int NProj = 10;
113 
114  String args=formatString("-i %s -o %s --sym %s --weight --thr %d -v 0",fnAngles.c_str(),fnVol.c_str(),fnSym.c_str(),Nthr);
115  String cmd=(String)"xmipp_reconstruct_fourier "+args;
116  if (system(cmd.c_str())==-1)
117  REPORT_ERROR(ERR_UNCLASSIFIED,"Cannot open shell");
118 
119  args=formatString("-i %s -o %s --mask circular %d -v 0",fnVol.c_str(),fnVol.c_str(),static_cast<int>(-xdim/2));
120  cmd=(String)"xmipp_transform_mask "+args;
121  if (system(cmd.c_str())==-1)
122  REPORT_ERROR(ERR_UNCLASSIFIED,"Cannot open shell");
123 
124  if (doProjMatch)
125  {
126 
127  for (int ind=0; ind<NProj; ind++)
128  {
129 
130  args=formatString("-i %s -o %s --sampling_rate %f --sym %s --method fourier 1 0.25 bspline --compute_neighbors --angular_distance -1 --experimental_images %s",fnVol.c_str(),fnGallery.c_str(),AngularSampling,fnSym.c_str(),fnAngles.c_str());
131  cmd=(String)"xmipp_angular_project_library "+args;
132  if (system(cmd.c_str())==-1)
133  REPORT_ERROR(ERR_UNCLASSIFIED,"Cannot open shell");
134 
135  args=formatString("-i %s -o %s --ref %s --Ri 0 --Ro %f --max_shift %f --append",fnAngles.c_str(),fnAngles.c_str(),fnGallery.c_str(),xdim/2,xdim/20);
136  cmd=(String)"xmipp_angular_projection_matching "+args;
137  if (system(cmd.c_str())==-1)
138  REPORT_ERROR(ERR_UNCLASSIFIED,"Cannot open shell");
139 
140  args=formatString("-i %s -o %s --sym %s --weight --thr %d -v 0",fnAngles.c_str(),fnVol.c_str(),fnSym.c_str(),Nthr);
141  cmd=(String)"xmipp_reconstruct_fourier "+args;
142  if (system(cmd.c_str())==-1)
143  REPORT_ERROR(ERR_UNCLASSIFIED,"Cannot open shell");
144 
145  args=formatString("-i %s -o %s --mask circular %d -v 0",fnVol.c_str(),fnVol.c_str(),static_cast<int>(-xdim/2));
146  cmd=(String)"xmipp_transform_mask "+args;
147  if (system(cmd.c_str())==-1)
148  REPORT_ERROR(ERR_UNCLASSIFIED,"Cannot open shell");
149  }
150  }
151 
152  args=formatString("-i %s -o %s --sym %s --weight --thr %d -v 0",fnAngles.c_str(),fnVol.c_str(),fnSym.c_str(),Nthr);
153  cmd=(String)"xmipp_reconstruct_fourier "+args;
154  if (system(cmd.c_str())==-1)
155  REPORT_ERROR(ERR_UNCLASSIFIED,"Cannot open shell");
156 
157  args=formatString("-i %s -o %s --mask circular %d -v 0",fnVol.c_str(),fnVol.c_str(),static_cast<int>(-xdim/2));
158  cmd=(String)"xmipp_transform_mask "+args;
159  if (system(cmd.c_str())==-1)
160  REPORT_ERROR(ERR_UNCLASSIFIED,"Cannot open shell");
161 
162  args=formatString("-i %s -o %s --select below 0 --substitute value 0 -v 0",fnVol.c_str(),fnVol.c_str());
163  cmd=(String)"xmipp_transform_threshold "+args;
164  if (system(cmd.c_str())==-1)
165  REPORT_ERROR(ERR_UNCLASSIFIED,"Cannot open shell");
166 }
ERR_UNCLASSIFIED
Just to locate unclassified errors.
Definition: xmipp_error.h:192
REPORT_ERROR
#define REPORT_ERROR(nerr, ErrormMsg)
Definition: xmipp_error.h:211
String
std::string String
Definition: xmipp_strings.h:34
formatString
String formatString(const char *format,...)
Definition: xmipp_strings.cpp:602

◆ run()

void ProgVolumeValidationPCA::run ( )
virtual

Run.

Reimplemented from XmippProgram.

Definition at line 223 of file volume_validate_pca.cpp.

224 {
225  show();
226  produceSideinfo();
227  std::stringstream ss;
228 
229  //Here we read the "mean" volume
230  //fnAngles=fnClasses;
231  //fnVol=fnClasses.removeAllExtensions();
232  //fnVol+=".vol";
233  //reconstruct();
234 
235  doProjMatch = true;
236 
237  for( int index = 0; index< NVols; index++)
238  {
239 
240  ss << index;
241  fnAngles=fnClasses.removeAllExtensions()+ss.str();
242  fnAngles+=".xmd";
243  ss.str(std::string());
244 
245  modifyAngles();
246 
247  ss << index;
248  fnVol=fnClasses.removeAllExtensions()+ss.str();
249  fnVol+=".vol";
250  ss.str(std::string());
251 
252  ss << index;
253  fnGallery=fnClasses.removeAllExtensions()+ss.str();
254  fnGallery+=".stk";
255  ss.str(std::string());
256 
257  //Reconstruct and ProjMatching of current
258  reconstruct();
259 
260  //Some clean up
261  FileName tempSampling = fnVol.removeAllExtensions();
262  tempSampling+="_sampling.xmd";
263  FileName tempDoc = fnAngles.removeAllExtensions();
264  tempDoc+=".doc";
265 
266  //String args=formatString(" %s %s %s %s",fnVol.c_str(),fnGallery.c_str(),tempSampling.c_str(),tempDoc.c_str());
267  //String cmd=(String)"rm "+args;
268  //system(cmd.c_str());
269  }
270 
271  evaluate();
272 }
FileName::removeAllExtensions
FileName removeAllExtensions() const
Definition: xmipp_filename.cpp:382
ProgVolumeValidationPCA::modifyAngles
void modifyAngles()
Modified angles.
Definition: volume_validate_pca.cpp:74
index
viol index
Definition: numerical_recipes.cpp:5197
ProgVolumeValidationPCA::produceSideinfo
void produceSideinfo()
Produce side info: fill arrays with relevant transformation matrices.
Definition: volume_validate_pca.cpp:66
ProgVolumeValidationPCA::reconstruct
void reconstruct()
Reconstruct volume with all projections.
Definition: volume_validate_pca.cpp:107
ProgVolumeValidationPCA::evaluate
void evaluate()
Modified angles.
Definition: volume_validate_pca.cpp:168

◆ show()

void ProgVolumeValidationPCA::show ( )

Show.

Definition at line 54 of file volume_validate_pca.cpp.

55 {
56  if (verbose > 0)
57  {
58  std::cout << "Input classes metadata : " << fnClasses << std::endl;
59  std::cout << "Number of intermediate volumes to generate : " << NVols << std::endl;
60  std::cout << "Number of classes to be used : " << NClasses << std::endl;
61  if (fnSym != "")
62  std::cout << "Symmetry for projections : " << fnSym << std::endl;
63  }
64 }
XmippProgram::verbose
int verbose
Verbosity level.
Definition: xmipp_program.h:161

Member Data Documentation

◆ doProjMatch

bool ProgVolumeValidationPCA::doProjMatch

Definition at line 55 of file volume_validate_pca.h.

◆ fnAngles

FileName ProgVolumeValidationPCA::fnAngles

Definition at line 47 of file volume_validate_pca.h.

◆ fnClasses

FileName ProgVolumeValidationPCA::fnClasses

Filenames

Definition at line 47 of file volume_validate_pca.h.

◆ fnSym

FileName ProgVolumeValidationPCA::fnSym

Definition at line 47 of file volume_validate_pca.h.

◆ mdClasses

MetaDataVec ProgVolumeValidationPCA::mdClasses

Definition at line 62 of file volume_validate_pca.h.

◆ mdReconstruction

MetaDataVec ProgVolumeValidationPCA::mdReconstruction

Definition at line 62 of file volume_validate_pca.h.

◆ NClasses

int ProgVolumeValidationPCA::NClasses

Number of classes to generate the intermediate volumes

Definition at line 53 of file volume_validate_pca.h.

◆ NVols

int ProgVolumeValidationPCA::NVols

Number of intermediate volumes to generate

Definition at line 50 of file volume_validate_pca.h.

The documentation for this class was generated from the following files:
Generated by   doxygen 1.8.13