Xmipp  v3.23.11-Nereus
Classes | Functions | Variables
Alignment of volumes with Normal modes
Collaboration diagram for Alignment of volumes with Normal modes:

Classes

class  ProgNmaAlignmentVol
 
class  ObjFunc_nma_alignment_vol
 

Functions

 ProgNmaAlignmentVol::ProgNmaAlignmentVol ()
 Empty constructor. More...
 
 ProgNmaAlignmentVol::~ProgNmaAlignmentVol ()
 Destructor. More...
 
void ProgNmaAlignmentVol::defineParams ()
 Define params. More...
 
void ProgNmaAlignmentVol::readParams ()
 Read arguments from command line. More...
 
void ProgNmaAlignmentVol::show ()
 Show. More...
 
FileName ProgNmaAlignmentVol::createDeformedPDB () const
 
double ProgNmaAlignmentVol::computeFitness (Matrix1D< double > &trial) const
 
bool ProgNmaAlignmentVol::updateBestFit (double fitness)
 
virtual void ProgNmaAlignmentVol::preProcess ()
 
virtual void ProgNmaAlignmentVol::processImage (const FileName &fnImg, const FileName &fnImgOut, const MDRow &rowIn, MDRow &rowOut)
 
virtual void ProgNmaAlignmentVol::finishProcessing ()
 
virtual void ProgNmaAlignmentVol::writeVolumeParameters (const FileName &fnImg)
 
virtual void ProgNmaAlignmentVol::createWorkFiles ()
 
 ObjFunc_nma_alignment_vol::ObjFunc_nma_alignment_vol (int _t, int _n=0)
 
 ObjFunc_nma_alignment_vol::~ObjFunc_nma_alignment_vol ()
 
double ObjFunc_nma_alignment_vol::eval (Vector v, int *nerror=nullptr)
 

Variables

bool ProgNmaAlignmentVol::MPIversion
 
bool ProgNmaAlignmentVol::resume
 
FileName ProgNmaAlignmentVol::fnPDB
 Reference atomic or pseudo-atomic structure in PDB format. More...
 
FileName ProgNmaAlignmentVol::fnOutPDB
 Output PDB. More...
 
FileName ProgNmaAlignmentVol::fnOutDir
 Output directory. More...
 
FileName ProgNmaAlignmentVol::fnModeList
 File with a list of mode filenames. More...
 
double ProgNmaAlignmentVol::sampling_rate
 Pixel size in Angstroms. More...
 
FileName ProgNmaAlignmentVol::fnmask
 Mask file for 2D masking of the projections of the deformed volume. More...
 
bool ProgNmaAlignmentVol::do_centerPDB
 Center the PDB structure. More...
 
bool ProgNmaAlignmentVol::do_FilterPDBVol
 Low-pass filter the volume from PDB. More...
 
double ProgNmaAlignmentVol::cutoff_LPfilter
 Low-pass cut-off frequency. More...
 
bool ProgNmaAlignmentVol::useFixedGaussian
 Use pseudo-atoms instead of atoms. More...
 
double ProgNmaAlignmentVol::sigmaGaussian
 Gaussian standard deviation for pseudo-atoms. More...
 
bool ProgNmaAlignmentVol::alignVolumes
 Align volumes. More...
 
double ProgNmaAlignmentVol::trustradius_scale
 Parameters required from the CONDOR optimization. More...
 
double ProgNmaAlignmentVol::rhoStartBase
 
double ProgNmaAlignmentVol::rhoEndBase
 
int ProgNmaAlignmentVol::niter
 
int ProgNmaAlignmentVol::tilt0
 
int ProgNmaAlignmentVol::tiltF
 
double ProgNmaAlignmentVol::frm_freq
 
int ProgNmaAlignmentVol::frm_shift
 
int ProgNmaAlignmentVol::rangen
 
Matrix1D< double > ProgNmaAlignmentVol::parameters
 
Matrix1D< double > ProgNmaAlignmentVol::trial
 
Matrix1D< double > ProgNmaAlignmentVol::trial_best
 
double ProgNmaAlignmentVol::fitness_min
 
int ProgNmaAlignmentVol::numberOfModes
 
int ProgNmaAlignmentVol::imgSize
 
FileName ProgNmaAlignmentVol::currentVolName
 
char ProgNmaAlignmentVol::nameTemplate [256]
 
ProgPdbConverterProgNmaAlignmentVol::progVolumeFromPDB
 
Image< double > ProgNmaAlignmentVol::V
 
Image< double > ProgNmaAlignmentVol::Vdeformed
 
MultidimArray< int > ProgNmaAlignmentVol::mask
 
FILE * ProgNmaAlignmentVol::AnglesShiftsAndScore
 
float ProgNmaAlignmentVol::Best_Angles_Shifts [6]
 
float ProgNmaAlignmentVol::fit_value
 
bool ProgNmaAlignmentVol::flip = false
 

Detailed Description

Function Documentation

◆ computeFitness()

double ProgNmaAlignmentVol::computeFitness ( Matrix1D< double > &  trial) const

Computes the fitness of a set of trial parameters

◆ createDeformedPDB()

FileName ProgNmaAlignmentVol::createDeformedPDB ( ) const

Create deformed PDB

Definition at line 164 of file nma_alignment_vol.cpp.

164  {
165  String program;
166  String arguments;
167  FileName fnRandom;
169  const char * randStr = fnRandom.c_str();
170 
171  deformPDB(fnPDB,formatString("%s_deformedPDB.pdb",randStr),fnModeList,trial);
172  program = "xmipp_volume_from_pdb";
173  arguments = formatString(
174  "-i %s_deformedPDB.pdb --size %i --sampling %f -v 0", randStr,
176 
177  if (do_centerPDB)
178  arguments.append(" --centerPDB ");
179 
180  if (useFixedGaussian) {
181  arguments.append(" --intensityColumn Bfactor --fixed_Gaussian ");
182  if (sigmaGaussian >= 0)
183  arguments += formatString("%f",sigmaGaussian);
184  }
185 
186  progVolumeFromPDB->read(arguments);
188 
189  if (do_FilterPDBVol) {
190  program = "xmipp_transform_filter";
191  arguments = formatString(
192  "-i %s_deformedPDB.vol --sampling %f --fourier low_pass %f -v 0",
193  randStr, sampling_rate, cutoff_LPfilter);
194  runSystem(program, arguments, false);
195  }
196 
197  return fnRandom;
198 }
virtual void read(int argc, const char **argv, bool reportErrors=true)
FileName fnOutDir
Output directory.
virtual int tryRun()
bool do_FilterPDBVol
Low-pass filter the volume from PDB.
double sampling_rate
Pixel size in Angstroms.
FileName fnPDB
Reference atomic or pseudo-atomic structure in PDB format.
FileName fnModeList
File with a list of mode filenames.
void initUniqueName(const char *templateStr="xmippTemp_XXXXXX", const String &fnDir="")
void deformPDB(const FileName &PDBin, const FileName &PDBout, const FileName &fnModeList, const Matrix1D< double > &trial)
void runSystem(const String &program, const String &arguments, bool useSystem)
bool do_centerPDB
Center the PDB structure.
bool useFixedGaussian
Use pseudo-atoms instead of atoms.
ProgPdbConverter * progVolumeFromPDB
double sigmaGaussian
Gaussian standard deviation for pseudo-atoms.
std::string String
Definition: xmipp_strings.h:34
Matrix1D< double > trial
String formatString(const char *format,...)
double cutoff_LPfilter
Low-pass cut-off frequency.

◆ createWorkFiles()

virtual void ProgNmaAlignmentVol::createWorkFiles ( )
inlineprotectedvirtual

Definition at line 183 of file nma_alignment_vol.h.

183  {
185  }
virtual void createWorkFiles(bool resume, MetaData *md)

◆ defineParams()

void ProgNmaAlignmentVol::defineParams ( )
virtual

Define params.

Reimplemented from XmippMetadataProgram.

Definition at line 46 of file nma_alignment_vol.cpp.

46  {
47  addUsageLine("Align volumes with an atomic or pseudo-atomic (from EM volume) structure by computing deformation amplitudes along normal modes");
48  addUsageLine("You may also align the two volumes during the process");
49  defaultComments["-i"].clear();
50  defaultComments["-i"].addComment("Metadata with volume filenames");
51  defaultComments["-o"].clear();
52  defaultComments["-o"].addComment("Metadata with output Euler angles, shifts, and deformation amplitudes");
54  addParamsLine(" --pdb <PDB_filename> : Reference atomic or pseudo-atomic structure in PDB format");
55  addParamsLine(" [--odir <outputDir=\".\">] : Output directory");
56  addParamsLine(" [--resume] : Resume processing");
57  addParamsLine(" [--opdb <PDB_filename=\"\">] : Output PDB which is the deformed input PDB");
58  addParamsLine("==Generation of the deformed reference volumes==");
59  addParamsLine(" --modes <filename> : File with a list of mode filenames");
60  addParamsLine(" [--sampling_rate <Ts=1>] : Pixel size in Angstroms");
61  addParamsLine(" [--filterVol <cutoff=15.>] : Filter the volume after deforming. If this option is used, the default cut-off is 15 A.");
62  addParamsLine(" [--centerPDB] : Center the PDB structure");
63  addParamsLine(" [--fixed_Gaussian <std=-1>] : For pseudo-atoms fixed_Gaussian must be used.");
64  addParamsLine(" : Default standard deviation <std> is read from the PDB file.");
65  addParamsLine(" [--trustradius_scale <s=1>] : Positive scaling factor to scale the initial trust region radius");
66  addParamsLine("==Combined elastic and rigid-body alignment==");
67  addParamsLine(" [--alignVolumes <frm_freq=0.25> <frm_shift=10>] : Align the deformed volume to the input volume before comparing, this is using frm method for volume alignment with frm_freq and frm_shift as parameters");
68  addParamsLine(" : You need to compile Xmipp with SHALIGNMENT support (see install.sh)");
69  addParamsLine(" [--mask <m=\"\">] : 3D masking of the projections of the deformed volume");
70  addParamsLine(" [--tilt_values <tilt0=-90> <tiltF=90>] : only use if you are trying to compensate for the missing wedge");
71  addParamsLine(" [--condor_params <rhoStartBase=250.> <rhoEndBase=50.> <niter=10000>] : parameters for the CONDOR optimiser (recommended to keep default)");
72  addParamsLine(" : rhoStartBase > 0 : the lower the better, yet the slower");
73  addParamsLine(" : rhoEndBase no specific rule, however it is better to keep it < 1000, if set very high we risk distortions");
74  addParamsLine(" : niter should be big enough to guarantee that the search converges to the right set of nma deformation amplitudes");
75  addExampleLine("xmipp_nma_alignment_vol -i volumes.xmd --pdb 2tbv.pdb --modes modelist.xmd --sampling_rate 3.2 -o output.xmd --resume");
76 }
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

◆ eval()

double ObjFunc_nma_alignment_vol::eval ( Vector  v,
int *  nerror = nullptr 
)
virtual

Implements ObjectiveFunction.

Definition at line 210 of file nma_alignment_vol.cpp.

210  {
212 
213  for (int i = 0; i < dim; i++) {
214  global_nma_vol_prog->trial(i) = X[i];
215  }
216 
218  const char * randStr = fnRandom.c_str();
219  double retval=std::numeric_limits<double>::max();
220 
221  FileName fnShiftsAngles = fnRandom + "_angles_shifts.txt";
222  const char * shifts_angles = fnShiftsAngles.c_str();
223 
225  String fnVolume2 = fnRandom + "_deformedPDB.vol";
226 
228  global_nma_vol_prog->flip = true;
229  runSystem("xmipp_transform_geometry",formatString("-i %s -o %s_currentvolume.vol --rotate_volume euler 0 90 0 -v 0",fnVolume1.c_str(),randStr));
230  fnVolume1 = fnRandom+"_currentvolume.vol";
231  }
232 
233  const char * Volume1 = fnVolume1.c_str();
234  const char * Volume2 = fnVolume2.c_str();
235 
236  int err;
237 
239  runSystem("xmipp_volume_align",formatString("--i1 %s --i2 %s --frm %f %d %d %d --store %s -v 0 ",
241  //first just see what is the score
242  global_nma_vol_prog->AnglesShiftsAndScore = fopen(shifts_angles, "r");
243  //fit_value is the 7th element in a single line CSV
244  for (int i = 0; i < 7; i++){
246  if (1!=err)
247  REPORT_ERROR(ERR_IO, "reading the fitness value was not successful");
248  }
250  retval = 1 + global_nma_vol_prog->fit_value;
251  }
252 
253  else{
254  global_nma_vol_prog->Vdeformed.read(formatString("%s_deformedPDB.vol",fnRandom.c_str()));
255  auto mask = (0 == XSIZE(global_nma_vol_prog->mask)) ? nullptr : &global_nma_vol_prog->mask;
257  //global_nma_vol_prog->V().printStats();
258  //global_nma_vol_prog->Vdeformed().printStats();
259  //std::cout << correlationIndex(global_nma_vol_prog->V(),global_nma_vol_prog->Vdeformed()) << std::endl;
260  }
261 
263  global_nma_vol_prog->AnglesShiftsAndScore = fopen(shifts_angles, "r");
264  for (int i = 0; i < 6; i++){
266  if (1!=err)
267  REPORT_ERROR(ERR_IO, "reading the angles and shifts was not successful");
268  }
270  }
271 
272  runSystem("rm", formatString("-rf %s* &", randStr));
273  //std::cout << global_nma_vol_prog->trial << " -> " << retval << std::endl;
274  return retval;
275 }
#define REPORT_ERROR(nerr, ErrormMsg)
Definition: xmipp_error.h:211
double correlationIndex(const MultidimArray< T > &x, const MultidimArray< T > &y, const MultidimArray< int > *mask=NULL, MultidimArray< double > *Contributions=NULL)
Input/Output general error.
Definition: xmipp_error.h:134
bool updateBestFit(double fitness)
void runSystem(const String &program, const String &arguments, bool useSystem)
#define i
FileName createDeformedPDB() const
#define XSIZE(v)
void max(Image< double > &op1, const Image< double > &op2)
bool alignVolumes
Align volumes.
std::string String
Definition: xmipp_strings.h:34
Matrix1D< double > trial
String formatString(const char *format,...)
int read(const FileName &name, DataMode datamode=DATA, size_t select_img=ALL_IMAGES, bool mapData=false, int mode=WRITE_READONLY)
ProgNmaAlignmentVol * global_nma_vol_prog
MultidimArray< int > mask
Image< double > Vdeformed

◆ finishProcessing()

void ProgNmaAlignmentVol::finishProcessing ( )
virtual

Write the final parameters.

Reimplemented from XmippMetadataProgram.

Definition at line 148 of file nma_alignment_vol.cpp.

148  {
150  rename(Rerunable::getFileName().c_str(), fn_out.c_str());
151 }
const FileName & getFileName() const

◆ ObjFunc_nma_alignment_vol()

ObjFunc_nma_alignment_vol::ObjFunc_nma_alignment_vol ( int  _t,
int  _n = 0 
)

Definition at line 277 of file nma_alignment_vol.cpp.

277  {
278 }

◆ preProcess()

void ProgNmaAlignmentVol::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 nma_alignment_vol.cpp.

130  {
132  SF.removeDisabled();
133  numberOfModes = SF.size();
134  // Get the size of the images in the selfile
135  imgSize = xdimOut;
136  // Set the pointer of the program to this object
137  global_nma_vol_prog = this;
138  //create some neededs files
139  createWorkFiles();
140  if (fnmask!="")
141  {
142  Image<double> aux;
143  aux.read(fnmask);
144  typeCast(aux(),mask);
145  }
146 }
FileName fnModeList
File with a list of mode filenames.
FileName fnmask
Mask file for 2D masking of the projections of the deformed volume.
void typeCast(const Matrix1D< T1 > &v1, Matrix1D< T2 > &v2)
Definition: matrix1d.h:1227
int read(const FileName &name, DataMode datamode=DATA, size_t select_img=ALL_IMAGES, bool mapData=false, int mode=WRITE_READONLY)
ProgNmaAlignmentVol * global_nma_vol_prog
MultidimArray< int > mask
virtual void createWorkFiles()

◆ processImage()

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

Assign NMA and Alignment parameters to a volume

Implements XmippMetadataProgram.

Definition at line 280 of file nma_alignment_vol.cpp.

281  {
282  static size_t imageCounter = 0;
283  ++imageCounter;
284 
285  ObjectiveFunction *of;
286 
287  int dim = numberOfModes;
288 
289  parameters.initZeros(dim);
290  V.read(fnImg);
291  currentVolName = fnImg;
292  sprintf(nameTemplate, "_node%d_img%lu_XXXXXX", rangen, (long unsigned int)imageCounter);
293 
294  trial.initZeros(dim);
295  trial_best.initZeros(dim);
296 
297  fitness_min=1000000.0;
298 
299  of = new ObjFunc_nma_alignment_vol(1, dim);
300 
301  of->xStart.setSize(dim);
302  for (int i = 0; i < dim; i++)
303  of->xStart[i] = 0.;
304 
305  double rhoStart=trustradius_scale*rhoStartBase;
306  double rhoEnd=trustradius_scale*rhoEndBase;
307 
308  CONDOR(rhoStart, rhoEnd, niter, of);
309 
311 
314 
315  writeVolumeParameters(fnImg);
316  if (fnOutPDB!="")
318  delete of;
319 }
FileName fnOutPDB
Output PDB.
#define VEC_XSIZE(m)
Definition: matrix1d.h:77
FileName fnPDB
Reference atomic or pseudo-atomic structure in PDB format.
FileName fnModeList
File with a list of mode filenames.
void deformPDB(const FileName &PDBin, const FileName &PDBout, const FileName &fnModeList, const Matrix1D< double > &trial)
#define i
void setSize(int _n)
Definition: Vector.cpp:112
Matrix1D< double > parameters
virtual void writeVolumeParameters(const FileName &fnImg)
double trustradius_scale
Parameters required from the CONDOR optimization.
void resize(size_t Xdim, bool copy=true)
Definition: matrix1d.h:410
Matrix1D< double > trial_best
void initZeros()
Definition: matrix1d.h:592
Matrix1D< double > trial
int read(const FileName &name, DataMode datamode=DATA, size_t select_img=ALL_IMAGES, bool mapData=false, int mode=WRITE_READONLY)
void CONDOR(double rhoStart, double rhoEnd, int niter, ObjectiveFunction *of, int nnode)
Definition: CNLSolver.cpp:75

◆ ProgNmaAlignmentVol()

ProgNmaAlignmentVol::ProgNmaAlignmentVol ( )

Empty constructor.

Definition at line 31 of file nma_alignment_vol.cpp.

31  : Rerunable("") {
32  rangen = 0;
33  resume = false;
34  currentVolName = "";
36  produces_an_output = true;
38  alignVolumes=false;
39 }
bool produces_an_output
Indicate that a unique final output is produced.
ProgPdbConverter * progVolumeFromPDB
bool alignVolumes
Align volumes.
bool each_image_produces_an_output
Indicate that an output is produced for each image in the input.
Rerunable(const FileName &fn)

◆ readParams()

void ProgNmaAlignmentVol::readParams ( )
virtual

Read arguments from command line.

Reimplemented from XmippMetadataProgram.

Definition at line 79 of file nma_alignment_vol.cpp.

79  {
81  fnPDB = getParam("--pdb");
82  fnOutPDB = getParam("--opdb");
83  fnOutDir = getParam("--odir");
84  Rerunable::setFileName(fnOutDir + "/nmaDone.xmd");
85  fnModeList = getParam("--modes");
86  resume = checkParam("--resume");
87  sampling_rate = getDoubleParam("--sampling_rate");
88  fnmask = getParam("--mask");
89  do_centerPDB = checkParam("--centerPDB");
90  do_FilterPDBVol = checkParam("--filterVol");
91  trustradius_scale = std::abs(getDoubleParam("--trustradius_scale"));
92  if (do_FilterPDBVol)
93  cutoff_LPfilter = getDoubleParam("--filterVol");
94  useFixedGaussian = checkParam("--fixed_Gaussian");
95  if (useFixedGaussian)
96  sigmaGaussian = getDoubleParam("--fixed_Gaussian");
97  alignVolumes=checkParam("--alignVolumes");
98  if (alignVolumes){
99  frm_freq = getDoubleParam("--alignVolumes",0);
100  frm_shift= getIntParam("--alignVolumes",1);
101  }
102  tilt0=getIntParam("--tilt_values",0);
103  tiltF=getIntParam("--tilt_values",1);
104 
105  rhoStartBase = getDoubleParam("--condor_params",0);
106  rhoEndBase = getDoubleParam("--condor_params",1);
107  niter = getIntParam("--condor_params",2);
108 
109 }
FileName fnOutPDB
Output PDB.
double getDoubleParam(const char *param, int arg=0)
FileName fnOutDir
Output directory.
bool do_FilterPDBVol
Low-pass filter the volume from PDB.
double sampling_rate
Pixel size in Angstroms.
FileName fnPDB
Reference atomic or pseudo-atomic structure in PDB format.
FileName fnModeList
File with a list of mode filenames.
void abs(Image< double > &op)
void setFileName(const FileName &fn)
FileName fnmask
Mask file for 2D masking of the projections of the deformed volume.
const char * getParam(const char *param, int arg=0)
bool do_centerPDB
Center the PDB structure.
double trustradius_scale
Parameters required from the CONDOR optimization.
bool useFixedGaussian
Use pseudo-atoms instead of atoms.
bool alignVolumes
Align volumes.
double sigmaGaussian
Gaussian standard deviation for pseudo-atoms.
double cutoff_LPfilter
Low-pass cut-off frequency.
bool checkParam(const char *param)
int getIntParam(const char *param, int arg=0)

◆ show()

void ProgNmaAlignmentVol::show ( )

Show.

Definition at line 111 of file nma_alignment_vol.cpp.

111  {
113  std::cout
114  << "Output directory: " << fnOutDir << std::endl
115  << "PDB: " << fnPDB << std::endl
116  << "Resume: " << resume << std::endl
117  << "Mode list: " << fnModeList << std::endl
118  << "Pixel size: " << sampling_rate << std::endl
119  << "Mask: " << fnmask << std::endl
120  << "Center PDB: " << do_centerPDB << std::endl
121  << "Filter PDB volume: " << do_FilterPDBVol << std::endl
122  << "Use pseudo-atoms: " << useFixedGaussian << std::endl
123  << "Pseudo-atom sigma: " << sigmaGaussian << std::endl
124  << "Trust-region scale: " << trustradius_scale << std::endl;
125 }
FileName fnOutDir
Output directory.
bool do_FilterPDBVol
Low-pass filter the volume from PDB.
double sampling_rate
Pixel size in Angstroms.
FileName fnPDB
Reference atomic or pseudo-atomic structure in PDB format.
FileName fnModeList
File with a list of mode filenames.
FileName fnmask
Mask file for 2D masking of the projections of the deformed volume.
bool do_centerPDB
Center the PDB structure.
double trustradius_scale
Parameters required from the CONDOR optimization.
bool useFixedGaussian
Use pseudo-atoms instead of atoms.
double sigmaGaussian
Gaussian standard deviation for pseudo-atoms.
void show() const override

◆ updateBestFit()

bool ProgNmaAlignmentVol::updateBestFit ( double  fitness)

Update the best fitness and the corresponding best trial

Definition at line 200 of file nma_alignment_vol.cpp.

200  {
201  if (fitness < fitness_min) {
203  trial_best = trial;
204  return true;
205  }
206  return false;
207 }
Matrix1D< double > trial_best
Matrix1D< double > trial
double fitness(double *p)

◆ writeVolumeParameters()

void ProgNmaAlignmentVol::writeVolumeParameters ( const FileName fnImg)
virtual

Write the parameters found for one image

Definition at line 321 of file nma_alignment_vol.cpp.

321  {
322  MetaDataVec md;
323  size_t objId = md.addObject();
324  md.setValue(MDL_IMAGE, fnImg, objId);
325  md.setValue(MDL_ENABLED, 1, objId);
326 
327  std::vector<double> vectortemp;
328  double energy=0;
329  for (int j = 0; j < numberOfModes; j++)
330  {
331  double lambdaj=parameters(j);
332  vectortemp.push_back(lambdaj);
333  energy+=lambdaj*lambdaj;
334  }
335 
336 
337  energy/=numberOfModes;
344  md.setValue(MDL_NMA, vectortemp, objId);
345  md.setValue(MDL_NMA_ENERGY, energy, objId);
346  md.setValue(MDL_MAXCC, 1-parameters(numberOfModes), objId);
348  md.setValue(MDL_ANGLE_Y, 90.0 , objId);
349  }
350  else{
351  md.setValue(MDL_ANGLE_Y, 0.0 , objId);
352  }
353 
355 }
Rotation angle of an image (double,degrees)
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.
Is this image enabled? (int [-1 or 1])
Matrix1D< double > parameters
const FileName & getFileName() const
bool setValue(const MDObject &mdValueIn, size_t id)
size_t addObject() override
Angle between y-axis and tilt-axis (double, degrees) for untilted micrographs.
Normal mode displacements (vector double)
Maximum cross-correlation for the image (double)
#define j
NMA energy contained in the NMA displacement vector.
void append(const FileName &outFile) const
Shift for the image in the Z axis (double)
ProgNmaAlignmentVol * global_nma_vol_prog
Shift for the image in the Y axis (double)
Name of an image (std::string)

◆ ~ObjFunc_nma_alignment_vol()

ObjFunc_nma_alignment_vol::~ObjFunc_nma_alignment_vol ( )
inline

Definition at line 203 of file nma_alignment_vol.h.

203 {};

◆ ~ProgNmaAlignmentVol()

ProgNmaAlignmentVol::~ProgNmaAlignmentVol ( )

Destructor.

Definition at line 41 of file nma_alignment_vol.cpp.

41  {
42  delete progVolumeFromPDB;
43 }
ProgPdbConverter * progVolumeFromPDB

Variable Documentation

◆ alignVolumes

bool ProgNmaAlignmentVol::alignVolumes

Align volumes.

Definition at line 82 of file nma_alignment_vol.h.

◆ AnglesShiftsAndScore

FILE* ProgNmaAlignmentVol::AnglesShiftsAndScore

Definition at line 137 of file nma_alignment_vol.h.

◆ Best_Angles_Shifts

float ProgNmaAlignmentVol::Best_Angles_Shifts[6]

Definition at line 138 of file nma_alignment_vol.h.

◆ currentVolName

FileName ProgNmaAlignmentVol::currentVolName

Definition at line 122 of file nma_alignment_vol.h.

◆ cutoff_LPfilter

double ProgNmaAlignmentVol::cutoff_LPfilter

Low-pass cut-off frequency.

Definition at line 73 of file nma_alignment_vol.h.

◆ do_centerPDB

bool ProgNmaAlignmentVol::do_centerPDB

Center the PDB structure.

Definition at line 67 of file nma_alignment_vol.h.

◆ do_FilterPDBVol

bool ProgNmaAlignmentVol::do_FilterPDBVol

Low-pass filter the volume from PDB.

Definition at line 70 of file nma_alignment_vol.h.

◆ fit_value

float ProgNmaAlignmentVol::fit_value

Definition at line 139 of file nma_alignment_vol.h.

◆ fitness_min

double ProgNmaAlignmentVol::fitness_min

Definition at line 113 of file nma_alignment_vol.h.

◆ flip

bool ProgNmaAlignmentVol::flip = false

Definition at line 142 of file nma_alignment_vol.h.

◆ fnmask

FileName ProgNmaAlignmentVol::fnmask

Mask file for 2D masking of the projections of the deformed volume.

Definition at line 64 of file nma_alignment_vol.h.

◆ fnModeList

FileName ProgNmaAlignmentVol::fnModeList

File with a list of mode filenames.

Definition at line 58 of file nma_alignment_vol.h.

◆ fnOutDir

FileName ProgNmaAlignmentVol::fnOutDir

Output directory.

Definition at line 55 of file nma_alignment_vol.h.

◆ fnOutPDB

FileName ProgNmaAlignmentVol::fnOutPDB

Output PDB.

Definition at line 52 of file nma_alignment_vol.h.

◆ fnPDB

FileName ProgNmaAlignmentVol::fnPDB

Reference atomic or pseudo-atomic structure in PDB format.

Definition at line 49 of file nma_alignment_vol.h.

◆ frm_freq

double ProgNmaAlignmentVol::frm_freq

Definition at line 94 of file nma_alignment_vol.h.

◆ frm_shift

int ProgNmaAlignmentVol::frm_shift

Definition at line 95 of file nma_alignment_vol.h.

◆ imgSize

int ProgNmaAlignmentVol::imgSize

Definition at line 119 of file nma_alignment_vol.h.

◆ mask

MultidimArray<int> ProgNmaAlignmentVol::mask

Definition at line 134 of file nma_alignment_vol.h.

◆ MPIversion

bool ProgNmaAlignmentVol::MPIversion

MPI version

Definition at line 43 of file nma_alignment_vol.h.

◆ nameTemplate

char ProgNmaAlignmentVol::nameTemplate[256]

Definition at line 125 of file nma_alignment_vol.h.

◆ niter

int ProgNmaAlignmentVol::niter

Definition at line 88 of file nma_alignment_vol.h.

◆ numberOfModes

int ProgNmaAlignmentVol::numberOfModes

Definition at line 116 of file nma_alignment_vol.h.

◆ parameters

Matrix1D<double> ProgNmaAlignmentVol::parameters

Definition at line 104 of file nma_alignment_vol.h.

◆ progVolumeFromPDB

ProgPdbConverter* ProgNmaAlignmentVol::progVolumeFromPDB

Definition at line 128 of file nma_alignment_vol.h.

◆ rangen

int ProgNmaAlignmentVol::rangen

Definition at line 101 of file nma_alignment_vol.h.

◆ resume

bool ProgNmaAlignmentVol::resume

Resume computations

Definition at line 46 of file nma_alignment_vol.h.

◆ rhoEndBase

double ProgNmaAlignmentVol::rhoEndBase

Definition at line 87 of file nma_alignment_vol.h.

◆ rhoStartBase

double ProgNmaAlignmentVol::rhoStartBase

Definition at line 86 of file nma_alignment_vol.h.

◆ sampling_rate

double ProgNmaAlignmentVol::sampling_rate

Pixel size in Angstroms.

Definition at line 61 of file nma_alignment_vol.h.

◆ sigmaGaussian

double ProgNmaAlignmentVol::sigmaGaussian

Gaussian standard deviation for pseudo-atoms.

Definition at line 79 of file nma_alignment_vol.h.

◆ tilt0

int ProgNmaAlignmentVol::tilt0

Definition at line 91 of file nma_alignment_vol.h.

◆ tiltF

int ProgNmaAlignmentVol::tiltF

Definition at line 91 of file nma_alignment_vol.h.

◆ trial

Matrix1D<double> ProgNmaAlignmentVol::trial

Definition at line 107 of file nma_alignment_vol.h.

◆ trial_best

Matrix1D<double> ProgNmaAlignmentVol::trial_best

Definition at line 110 of file nma_alignment_vol.h.

◆ trustradius_scale

double ProgNmaAlignmentVol::trustradius_scale

Parameters required from the CONDOR optimization.

Definition at line 85 of file nma_alignment_vol.h.

◆ useFixedGaussian

bool ProgNmaAlignmentVol::useFixedGaussian

Use pseudo-atoms instead of atoms.

Definition at line 76 of file nma_alignment_vol.h.

◆ V

Image<double> ProgNmaAlignmentVol::V

Definition at line 131 of file nma_alignment_vol.h.

◆ Vdeformed

Image<double> ProgNmaAlignmentVol::Vdeformed

Definition at line 131 of file nma_alignment_vol.h.