ProSHADE  0.7.5.4 (MAR 2021)
Protein Shape Detection
ProSHADE_settings Class Reference

This class stores all the settings and is passed to the executive classes instead of a multitude of parameters. More...

#include <ProSHADE_settings.hpp>

Public Member Functions

void determineBandwidthFromAngle (proshade_double uncertainty)
 This function determines the bandwidth for the spherical harmonics computation from the allowed rotation function angle uncertainty. More...
 
void determineBandwidth (proshade_unsign circumference)
 This function determines the bandwidth for the spherical harmonics computation. More...
 
void determineSphereDistances (proshade_single maxMapRange)
 This function determines the sphere distances for sphere mapping. More...
 
void determineIntegrationOrder (proshade_single maxMapRange)
 This function determines the integration order for the between spheres integration. More...
 
void determineAllSHValues (proshade_unsign xDim, proshade_unsign yDim, proshade_unsign zDim, proshade_double xDimAngs, proshade_double yDimAngs, proshade_double zDimAngs)
 This function determines all the required values for spherical harmonics computation. More...
 
void setVariablesLeftOnAuto (void)
 Function to determine general values that the user left on auto-determination.
 
 ProSHADE_settings (void)
 Contructor for the ProSHADE_settings class. More...
 
 ProSHADE_settings (ProSHADE_Task task)
 Contructor for the ProSHADE_settings class for particular task. More...
 
 ~ProSHADE_settings (void)
 Destructor for the ProSHADE_settings class. More...
 
void addStructure (std::string structure)
 Adds a structure file name to the appropriate variable. More...
 
void setResolution (proshade_single resolution)
 Sets the requested resolution in the appropriate variable. More...
 
void setPDBBFactor (proshade_double newBF)
 Sets the requested B-factor value for PDB files in the appropriate variable. More...
 
void setNormalisation (bool normalise)
 Sets the requested map normalisation value in the appropriate variable. More...
 
void setMapInversion (bool mInv)
 Sets the requested map inversion value in the appropriate variable. More...
 
void setVerbosity (proshade_signed verbosity)
 Sets the requested verbosity in the appropriate variable. More...
 
void setMaskBlurFactor (proshade_single blurFac)
 Sets the requested map blurring factor in the appropriate variable. More...
 
void setMaskIQR (proshade_single noIQRs)
 Sets the requested number of IQRs for masking threshold in the appropriate variable. More...
 
void setMasking (bool mask)
 Sets the requested map masking decision value in the appropriate variable. More...
 
void setCorrelationMasking (bool corMask)
 Sets the requested map masking type in the appropriate variable. More...
 
void setTypicalNoiseSize (proshade_single typNoi)
 Sets the requested "fake" half-map kernel in the appropriate variable. More...
 
void setMinimumMaskSize (proshade_single minMS)
 Sets the requested minimum mask size. More...
 
void setMaskSaving (bool savMsk)
 Sets whether the mask should be saved. More...
 
void setMaskFilename (std::string mskFln)
 Sets where the mask should be saved. More...
 
void setMapReboxing (bool reBx)
 Sets whether re-boxing needs to be done in the appropriate variable. More...
 
void setBoundsSpace (proshade_single boundsExSp)
 Sets the requested number of angstroms for extra space in re-boxing in the appropriate variable. More...
 
void setBoundsThreshold (proshade_signed boundsThres)
 Sets the threshold for number of indices difference acceptable to make index sizes same in the appropriate variable. More...
 
void setSameBoundaries (bool sameB)
 Sets whether same boundaries should be used in the appropriate variable. More...
 
void setOutputFilename (std::string oFileName)
 Sets the requested output file name in the appropriate variable. More...
 
void setMapResolutionChange (bool mrChange)
 Sets the requested map resolution change decision in the appropriate variable. More...
 
void setMapResolutionChangeTriLinear (bool mrChange)
 Sets the requested map resolution change decision using tri-linear interpolation in the appropriate variable. More...
 
void setMapCentering (bool com)
 Sets the requested map centering decision value in the appropriate variable. More...
 
void setExtraSpace (proshade_single exSpace)
 Sets the requested map extra space value in the appropriate variable. More...
 
void setBandwidth (proshade_unsign band)
 Sets the requested spherical harmonics bandwidth in the appropriate variable. More...
 
void setSphereDistances (proshade_single sphDist)
 Sets the requested distance between spheres in the appropriate variable. More...
 
void setIntegrationOrder (proshade_unsign intOrd)
 Sets the requested order for the Gauss-Legendre integration in the appropriate variable. More...
 
void setTaylorSeriesCap (proshade_unsign tayCap)
 Sets the requested Taylor series cap for the Gauss-Legendre integration in the appropriate variable. More...
 
void setProgressiveSphereMapping (bool progSphMap)
 Sets the requested sphere mapping value settings approach in the appropriate variable. More...
 
void setEnergyLevelsComputation (bool enLevDesc)
 Sets whether the energy level distance descriptor should be computed. More...
 
void setTraceSigmaComputation (bool trSigVal)
 Sets whether the trace sigma distance descriptor should be computed. More...
 
void setRotationFunctionComputation (bool rotfVal)
 Sets whether the rotation function distance descriptor should be computed. More...
 
void setPeakNeighboursNumber (proshade_unsign pkS)
 Sets the number of neighbour values that have to be smaller for an index to be considered a peak. More...
 
void setPeakNaiveNoIQR (proshade_double noIQRs)
 Sets the number of IQRs from the median for threshold height a peak needs to be considered a peak. More...
 
void setPhaseUsage (bool phaseUsage)
 Sets whether the phase information will be used. More...
 
void setEnLevShellWeight (proshade_double mPower)
 Sets the weight of shell position for the energy levels computation. More...
 
void setGroupingSmoothingFactor (proshade_double smFact)
 Sets the grouping smoothing factor into the proper variable. More...
 
void setMissingPeakThreshold (proshade_double mpThres)
 Sets the threshold for starting the missing peaks procedure. More...
 
void setAxisComparisonThreshold (proshade_double axThres)
 Sets the threshold for matching symmetry axes. More...
 
void setAxisComparisonThresholdBehaviour (bool behav)
 Sets the automatic symmetry axis tolerance decreasing. More...
 
void setMinimumPeakForAxis (proshade_double minSP)
 Sets the minimum peak height for symmetry axis to be considered. More...
 
void setRecommendedSymmetry (std::string val)
 Sets the ProSHADE detected symmetry type. More...
 
void setRecommendedFold (proshade_unsign val)
 Sets the ProSHADE detected symmetry fold. More...
 
void setRequestedSymmetry (std::string val)
 Sets the user requested symmetry type. More...
 
void setRequestedFold (proshade_unsign val)
 Sets the user requested symmetry fold. More...
 
void setDetectedSymmetry (proshade_double *sym)
 Sets the final detected symmetry axes information. More...
 
void setOverlaySaveFile (std::string filename)
 Sets the filename to which the overlay structure is to be save into. More...
 
void setOverlayJsonFile (std::string filename)
 Sets the filename to which the overlay operations are to be save into. More...
 
void setSymmetryRotFunPeaks (bool rotFunPeaks)
 Sets the symmetry detection algorithm type. More...
 
void setBicubicInterpolationSearch (bool bicubPeaks)
 Sets the bicubic interpolation on peaks. More...
 
void setMaxSymmetryFold (proshade_unsign maxFold)
 Sets the maximum symmetry fold (well, the maximum prime symmetry fold). More...
 
void getCommandLineParams (int argc, char **argv)
 This function parses the command line arguments into the settings object. More...
 
void printSettings (void)
 This function prints the current values in the settings object. More...
 

Public Attributes

ProSHADE_Task task
 This custom type variable determines which task to perfom (i.e. symmetry detection, distances computation, etc.).
 
std::vector< std::string > inputFiles
 This vector contains the filenames of all input structure files.
 
bool forceP1
 Should the P1 spacegroup be forced on the input PDB files?
 
bool removeWaters
 Should all waters be removed from input PDB files?
 
bool firstModelOnly
 Shoud only the first PDB model be used, or should all models be used?
 
proshade_single requestedResolution
 The resolution to which the calculations are to be done.
 
bool changeMapResolution
 Should maps be re-sampled to obtain the required resolution?
 
bool changeMapResolutionTriLinear
 Should maps be re-sampled to obtain the required resolution?
 
proshade_double pdbBFactorNewVal
 Change all PDB B-factors to this value (for smooth maps).
 
proshade_unsign maxBandwidth
 The bandwidth of spherical harmonics decomposition for the largest sphere.
 
proshade_double rotationUncertainty
 Alternative to bandwidth - the angle in degrees to which the rotation function accuracy should be computed.
 
bool usePhase
 If true, the full data will be used, if false, Patterson maps will be used instead and phased data will be converted to them. Also, only half of the spherical harmonics bands will be necessary as odd bands have to be 0 for Patterson maps.
 
proshade_single maxSphereDists
 The distance between spheres in spherical mapping for the largest sphere.
 
proshade_unsign integOrder
 The order required for full Gauss-Legendre integration between the spheres.
 
proshade_unsign taylorSeriesCap
 The max limit on the Taylor series expansion done for the abscissas of the Gauss-Legendre integration.
 
bool normaliseMap
 Should the map be normalised to mean 0 sd 1?
 
bool invertMap
 Should the map be inverted? Only use this if you think you have the wrong hand in your map.
 
proshade_single blurFactor
 This is the amount by which B-factors should be increased to create the blurred map for masking.
 
proshade_single maskingThresholdIQRs
 Number of inter-quartile ranges from the median to be used for thresholding the blurred map for masking.
 
bool maskMap
 Should the map be masked from noise?
 
bool useCorrelationMasking
 Should the blurring masking (false) or the correlation masking (true) be used?
 
proshade_single halfMapKernel
 This value in Angstrom will be used as the kernel for the "fake half-map" computation.
 
proshade_single correlationKernel
 This value in Angstrom will be used as the kernel for the map-FHM correlation computation.
 
bool saveMask
 Should the mask be saved?
 
std::string maskFileName
 The filename to which mask should be saved.
 
bool reBoxMap
 This switch decides whether re-boxing is needed.
 
proshade_single boundsExtraSpace
 The number of extra angstroms to be added to all re-boxing bounds just for safety.
 
proshade_signed boundsSimilarityThreshold
 Number of indices which can be added just to make sure same size in indices is achieved.
 
bool useSameBounds
 Switch to say that the same boundaries as used for the first should be used for all input maps.
 
proshade_signed * forceBounds
 These will be the boundaries to be forced upon the map.
 
bool moveToCOM
 Logical value stating whether the structure should be moved to have its Centre Of Mass (COM) in the middle.
 
proshade_single addExtraSpace
 If this value is non-zero, this many angstroms of empty space will be added to the internal map.
 
bool progressiveSphereMapping
 If true, each shell will have its own angular resolution dependent on the actual number of map points which are available to it. If false, all shells will have the same settings.
 
std::string outName
 The file name where the output structure(s) should be saved.
 
bool computeEnergyLevelsDesc
 If true, the energy levels descriptor will be computed, otherwise all its computations will be omitted.
 
proshade_double enLevMatrixPowerWeight
 If RRP matrices shell position is to be weighted by putting the position as an exponent, this variable sets the exponent. Set to 0 for no weighting.
 
bool computeTraceSigmaDesc
 If true, the trace sigma descriptor will be computed, otherwise all its computations will be omitted.
 
bool computeRotationFuncDesc
 If true, the rotation function descriptor will be computed, otherwise all its computations will be omitted.
 
proshade_unsign peakNeighbours
 Number of points in any direction that have to be lower than the considered index in order to consider this index a peak.
 
proshade_double noIQRsFromMedianNaivePeak
 When doing peak searching, how many IQRs from the median the threshold for peak height should be (in terms of median of non-peak values).
 
proshade_double smoothingFactor
 This factor decides how small the group sizes should be - larger factor means more smaller groups.
 
proshade_double symMissPeakThres
 Percentage of peaks that could be missing that would warrant starting the missing peaks search procedure.
 
proshade_double axisErrTolerance
 Allowed error on vector axis in in dot product ( acos ( 1 - axErr ) is the allowed difference in radians ).
 
bool axisErrToleranceDefault
 
proshade_double minSymPeak
 Minimum average peak for symmetry axis to be considered as "real".
 
std::string recommendedSymmetryType
 The symmetry type that ProSHADE finds the best fitting for the structure. Possible values are "" for none, "C" for cyclic, "D" for Dihedral, "T" for Tetrahedral, "O" for Octahedral and "I" for Icosahedral. C and D types also have fold value associated.
 
proshade_unsign recommendedSymmetryFold
 The fold of the recommended symmetry C or D type, 0 otherwise.
 
std::string requestedSymmetryType
 The symmetry type requested by the user. Allowed values are C, D, T, O and I.
 
proshade_unsign requestedSymmetryFold
 The fold of the requested symmetry (only applicable to C and D symmetry types).
 
bool usePeakSearchInRotationFunctionSpace
 This variable switch decides whether symmetry detection will be done using peak search in rotation function or using the angle-axis sperical space.
 
bool useBiCubicInterpolationOnPeaks
 This variable switch decides whether best symmetry is detected from peak indices, or whether bicubic interpolation is done to seatch for better axis between indices.
 
proshade_unsign maxSymmetryFold
 The highest symmetry fold to search for.
 
std::string overlayStructureName
 The filename to which the rotated and translated moving structure is to be saved.
 
std::string rotTrsJSONFile
 The filename to which the rotation and translation operations are to be saved into.
 
proshade_signed verbose
 Should the software report on the progress, or just be quiet? Value between -1 (nothing) and 4 (loud)
 
std::vector< proshade_double * > detectedSymmetry
 The vector of detected symmetry axes.
 
std::vector< std::vector< proshade_double > > allDetectedCAxes
 The vector of all detected cyclic symmetry axes.
 
std::vector< std::vector< proshade_unsign > > allDetectedDAxes
 The vector of all detected dihedral symmetry axes indices in allDetectedCAxes.
 
std::vector< proshade_unsign > allDetectedTAxes
 The vector of all detected tetrahedral symmetry axes indices in allDetectedCAxes.
 
std::vector< proshade_unsign > allDetectedOAxes
 The vector of all detected octahedral symmetry axes indices in allDetectedCAxes.
 
std::vector< proshade_unsign > allDetectedIAxes
 The vector of all detected icosahedral symmetry axes indices in allDetectedCAxes.
 

Detailed Description

This class stores all the settings and is passed to the executive classes instead of a multitude of parameters.

The ProSHADE_settings class is a simple way of keeping all the settings together and easy to set by the user. Its constructor sets it to the default settings, so that if the user does not want to change these, he just needs to pass the object to the executing class and all is done.

Definition at line 88 of file ProSHADE_settings.hpp.

Constructor & Destructor Documentation

◆ ProSHADE_settings() [1/2]

ProSHADE_settings::ProSHADE_settings ( void  )

Contructor for the ProSHADE_settings class.

This is the generic constructor used in cases where the settings object will be filled based on run-time determined values. If you know the specific task to be done, it is recommended to use the constructor which takes the task as argument, so that the default values are set specifically for the task at hand.

Definition at line 34 of file ProSHADE.cpp.

35 {
36  //================================================ Settings regarding the task at hand
37  this->task = NA;
38 
39  //================================================ Settings regarding input files
40  this->forceP1 = true;
41  this->removeWaters = true;
42  this->firstModelOnly = true;
43 
44  //================================================ Settings regarding the resolution of calculations
45  this->requestedResolution = -1.0;
46  this->changeMapResolution = false;
47  this->changeMapResolutionTriLinear = false;
48 
49  //================================================ Settings regarding the PDB B-factor change
50  this->pdbBFactorNewVal = -1.0;
51 
52  //================================================ Settings regarding the bandwidth of calculations
53  this->maxBandwidth = 0;
54  this->rotationUncertainty = 0;
55 
56  //================================================ Settings regarding the phase
57  this->usePhase = true;
58 
59  //================================================ Settings regarding the spheres
60  this->maxSphereDists = 0.0;
61 
62  //================================================ Settings regarding the Gauss-Legendre integration
63  this->integOrder = 0;
64  this->taylorSeriesCap = 10;
65 
66  //================================================ Settings regarding map normalisation
67  this->normaliseMap = false;
68 
69  //================================================ Settings regarding map inversion
70  this->invertMap = false;
71 
72  //================================================ Settings regarding map masking
73  this->blurFactor = 350.0;
74  this->maskingThresholdIQRs = 3.0;
75  this->maskMap = false;
76  this->useCorrelationMasking = false;
77  this->halfMapKernel = 0.0;
78  this->correlationKernel = 0.0;
79  this->saveMask = false;
80  this->maskFileName = "maskFile";
81 
82  //================================================ Settings regarding re-boxing
83  this->reBoxMap = false;
84  this->boundsExtraSpace = 3.0;
85  this->boundsSimilarityThreshold = 0;
86  this->useSameBounds = false;
87  this->forceBounds = new proshade_signed [6];
88 
89  //================================================ Settings regarding COM
90  this->moveToCOM = false;
91 
92  //================================================ Settings regarding extra cell space
93  this->addExtraSpace = 10.0;
94 
95  //================================================ Settings regarding shell settings
96  this->progressiveSphereMapping = false;
97 
98  //================================================ Settings regarding output file name
99  this->outName = "reBoxed";
100 
101  //================================================ Settings regarding distances computation
102  this->computeEnergyLevelsDesc = true;
103  this->computeTraceSigmaDesc = true;
104  this->computeRotationFuncDesc = true;
105  this->enLevMatrixPowerWeight = 1.0;
106 
107  //================================================ Settings regarding peak searching
108  this->peakNeighbours = 1;
109  this->noIQRsFromMedianNaivePeak = -999.9;
110 
111  //================================================ Settings regarding 1D grouping
112  this->smoothingFactor = 15.0;
113 
114  //================================================ Settings regarding the symmetry detection
116  this->useBiCubicInterpolationOnPeaks = true;
117  this->maxSymmetryFold = 30;
118  this->symMissPeakThres = 0.3;
119  this->axisErrTolerance = 0.01;
120  this->axisErrToleranceDefault = true;
121  this->minSymPeak = 0.3;
122  this->recommendedSymmetryType = "";
123  this->recommendedSymmetryFold = 0;
124  this->requestedSymmetryType = "";
125  this->requestedSymmetryFold = 0;
126  this->detectedSymmetry.clear ( );
127 
128  //================================================ Settings regarding the structure overlay
129  this->overlayStructureName = "movedStructure";
130  this->rotTrsJSONFile = "movedStructureOperations.json";
131 
132  //================================================ Settings regarding verbosity of the program
133  this->verbose = 1;
134 
135  //================================================ Done
136 
137 }

◆ ProSHADE_settings() [2/2]

ProSHADE_settings::ProSHADE_settings ( ProSHADE_Task  taskToPerform)

Contructor for the ProSHADE_settings class for particular task.

This is the generic constructor used in cases where the settings object will be filled based on run-time determined values. If you know the specific task to be done, it is recommended to use the constructor which takes the task as argument, so that the default values are set specifically for the task at hand.

Parameters
[in]taskToPerformThe task that should be performed by ProSHADE.

Definition at line 147 of file ProSHADE.cpp.

148 {
149  //================================================ Settings regarding the task at hand
150  this->task = taskToPerform;
151 
152  //================================================ Settings regarding input files
153  this->forceP1 = true;
154  this->removeWaters = true;
155  this->firstModelOnly = true;
156 
157  //================================================ Settings regarding the resolution of calculations
158  this->requestedResolution = -1.0;
159  this->changeMapResolution = false;
160  this->changeMapResolutionTriLinear = false;
161 
162  //================================================ Settings regarding the PDB B-factor change
163  this->pdbBFactorNewVal = -1.0;
164 
165  //================================================ Settings regarding the bandwidth of calculations
166  this->maxBandwidth = 0;
167  this->rotationUncertainty = 0;
168 
169  //================================================ Settings regarding the phase
170  this->usePhase = true;
171 
172  //================================================ Settings regarding the spheres
173  this->maxSphereDists = 0.0;
174 
175  //================================================ Settings regarding the Gauss-Legendre integration
176  this->integOrder = 0;
177  this->taylorSeriesCap = 10;
178 
179  //================================================ Settings regarding map normalisation
180  this->normaliseMap = false;
181 
182  //================================================ Settings regarding map inversion
183  this->invertMap = false;
184 
185  //================================================ Settings regarding map masking
186  this->blurFactor = 350.0;
187  this->maskingThresholdIQRs = 3.0;
188  this->maskMap = false;
189  this->useCorrelationMasking = false;
190  this->halfMapKernel = 0.0;
191  this->correlationKernel = 0.0;
192  this->saveMask = false;
193  this->maskFileName = "maskFile";
194  this->detectedSymmetry.clear ( );
195 
196  //================================================ Settings regarding re-boxing
197  this->reBoxMap = false;
198  this->boundsExtraSpace = 3.0;
199  this->boundsSimilarityThreshold = 0;
200  this->useSameBounds = false;
201  this->forceBounds = new proshade_signed [6];
202 
203  //================================================ Settings regarding extra cell space
204  this->addExtraSpace = 10.0;
205 
206  //================================================ Settings regarding shell settings
207  this->progressiveSphereMapping = false;
208 
209  //================================================ Settings regarding output file name
210  this->outName = "reBoxed";
211 
212  //================================================ Settings regarding distances computation
213  this->computeEnergyLevelsDesc = true;
214  this->computeTraceSigmaDesc = true;
215  this->computeRotationFuncDesc = true;
216  this->enLevMatrixPowerWeight = 1.0;
217 
218  //================================================ Settings regarding peak searching
219  this->peakNeighbours = 1;
220  this->noIQRsFromMedianNaivePeak = -999.9;
221 
222  //================================================ Settings regarding 1D grouping
223  this->smoothingFactor = 15.0;
224 
225  //================================================ Settings regarding the symmetry detection
227  this->useBiCubicInterpolationOnPeaks = true;
228  this->maxSymmetryFold = 30;
229  this->symMissPeakThres = 0.3;
230  this->axisErrTolerance = 0.01;
231  this->axisErrToleranceDefault = true;
232  this->minSymPeak = 0.3;
233  this->recommendedSymmetryType = "";
234  this->recommendedSymmetryFold = 0;
235  this->requestedSymmetryType = "";
236  this->requestedSymmetryFold = 0;
237 
238  //================================================ Settings regarding the structure overlay
239  this->overlayStructureName = "movedStructure";
240  this->rotTrsJSONFile = "movedStructureOperations.json";
241 
242  //================================================ Settings regarding verbosity of the program
243  this->verbose = 1;
244 
245  //================================================ Task specific settings
246  switch ( this->task )
247  {
248  case NA:
249  std::cerr << std::endl << "=====================" << std::endl << "!! ProSHADE ERROR !!" << std::endl << "=====================" << std::endl << std::flush;
250  std::cerr << "Error Code : " << "E000014" << std::endl << std::flush;
251  std::cerr << "ProSHADE version : " << __PROSHADE_VERSION__ << std::endl << std::flush;
252  std::cerr << "File : " << "ProSHADE.cpp" << std::endl << std::flush;
253  std::cerr << "Line : " << 97 << std::endl << std::flush;
254  std::cerr << "Function : " << "ProSHADE_settings (Task) constructor" << std::endl << std::flush;
255  std::cerr << "Message : " << "No task has been specified for task specific constructor." << std::endl << std::flush;
256  std::cerr << "Further information : " << "This ProSHADE_settings class constructor is intended to\n : set the internal variables to default value given a\n : particular taks. By supplying this task as NA, this beats\n : the purpose of the constructor. Please use the\n : non-argumental constructor if task is not yet known." << std::endl << std::endl << std::flush;
258  exit ( EXIT_FAILURE );
259  break;
260 
261  case Symmetry:
262  this->requestedResolution = 6.0;
263  this->pdbBFactorNewVal = 80.0;
264  this->changeMapResolution = true;
265  this->maskMap = false;
266  this->moveToCOM = true;
267  this->normaliseMap = true;
268  this->reBoxMap = false;
269  break;
270 
271  case Distances:
272  this->changeMapResolution = false;
273  this->maskMap = false;
274  this->moveToCOM = true;
275  this->reBoxMap = false;
276  break;
277 
278  case OverlayMap:
279  this->requestedResolution = 8.0;
280  this->changeMapResolution = true;
281  this->maskMap = false;
282  this->moveToCOM = false;
283  this->normaliseMap = false;
284  this->reBoxMap = false;
285  break;
286 
287  case MapManip:
288  this->changeMapResolution = false;
289  this->maskMap = true;
290  this->moveToCOM = false;
291  break;
292  }
293 
294  //================================================ Done
295 
296 }

◆ ~ProSHADE_settings()

ProSHADE_settings::~ProSHADE_settings ( void  )

Destructor for the ProSHADE_settings class.

This destructor is responsible for releasing all memory used by the settings object

Definition at line 302 of file ProSHADE.cpp.

303 {
304  //================================================ Release boundaries variable
305  delete[] this->forceBounds;
306 
307  //================================================ Release symmetry axes
308  if ( this->detectedSymmetry.size() > 0 ) { for ( proshade_unsign it = 0; it < static_cast<proshade_unsign> ( this->detectedSymmetry.size() ); it++ ) { if ( this->detectedSymmetry.at(it) != NULL ) { delete[] this->detectedSymmetry.at(it); } } }
309 
310  //================================================ Done
311 
312 }

Member Function Documentation

◆ addStructure()

void ProSHADE_settings::addStructure ( std::string  structure)

Adds a structure file name to the appropriate variable.

This function takes a string defining the filename of a structure to be processed and adds it to the list of structures to be processed.

Parameters
[in]structureString file name to be added to the list of structures to process.

Definition at line 341 of file ProSHADE.cpp.

342 {
343  //================================================ Use C++ version independent vector processing
344  ProSHADE_internal_misc::addToStringVector ( &( this->inputFiles ), structure );
345 
346  //================================================ Done
347  return ;
348 
349 }

◆ determineAllSHValues()

void ProSHADE_settings::determineAllSHValues ( proshade_unsign  xDim,
proshade_unsign  yDim,
proshade_unsign  zDim,
proshade_double  xDimAngs,
proshade_double  yDimAngs,
proshade_double  zDimAngs 
)

This function determines all the required values for spherical harmonics computation.

This function takes the maximum dimension size (in indices) and uses the settings pre-set by the user to set up the sphherical harmonics bandwidth, sphere sampling, sphere placement and spacing as well as the Gauss-Legendre integration order. This is either done using the user set values (if given), or using automated algorithm which only requires the resolution and max dimension.

Note that this function will use the resolution value to modify the values to be appropriate for the resolution supplied and not necessarily for the map sampling given.

Parameters
[in]xDimThe size of the x axis dimension in indices.
[in]yDimThe size of the y axis dimension in indices.
[in]zDimThe size of the z axis dimension in indices.
[in]xDimAngsThe size of the x-axis in Angstroms.
[in]yDimAngsThe size of the y-axis in Angstroms.
[in]zDimAngsThe size of the z-axis in Angstroms.
Warning
Because the automated algorithm decides the values based on the first structure size, by using it one gives up on the idea that DIST(A,B) == DIST(B,A). If this is important, then the user should set all of these values manually to the settings object to avoid this issue.

Definition at line 1320 of file ProSHADE.cpp.

1321 {
1322  //================================================ Print progress message
1323  ProSHADE_internal_messages::printProgressMessage ( this->verbose, 1, "Preparing spherical harmonics environment." );
1324 
1325  //================================================ Modify dims by resolution
1326  proshade_unsign theoXDim = std::ceil ( xDimAngs / ( this->requestedResolution / 2.0 ) );
1327  proshade_unsign theoYDim = std::ceil ( yDimAngs / ( this->requestedResolution / 2.0 ) );
1328  proshade_unsign theoZDim = std::ceil ( zDimAngs / ( this->requestedResolution / 2.0 ) );
1329 
1330  //================================================ Find maximum circumference
1331  proshade_unsign maxDim = std::max ( theoXDim, std::max ( theoYDim, theoZDim ) );
1332  proshade_unsign minDim = std::min ( theoXDim, std::min ( theoYDim, theoZDim ) );
1333  proshade_unsign midDim = 0;
1334  if ( ( xDim < maxDim ) && ( xDim > minDim ) ) { midDim = theoXDim; }
1335  else if ( ( yDim < maxDim ) && ( yDim > minDim ) ) { midDim = theoYDim; }
1336  else { midDim = theoZDim; }
1337 
1338  proshade_unsign circ = ( maxDim ) + ( midDim );
1339 
1340  //================================================ Bandwidth
1341  if ( this->rotationUncertainty > 0.0 ) { this->determineBandwidthFromAngle ( this->rotationUncertainty ); }
1342  else { this->determineBandwidth ( circ ); }
1343 
1344  //================================================ Find maximum diagonal in Angstroms
1345  proshade_single maxDiag = std::sqrt ( std::pow ( static_cast<proshade_single> ( maxDim ) * ( this->requestedResolution / 2.0 ), 2.0 ) +
1346  std::pow ( static_cast<proshade_single> ( midDim ) * ( this->requestedResolution / 2.0 ), 2.0 ) );
1347 
1348  //================================================ Sphere distances
1349  this->determineSphereDistances ( maxDiag );
1350 
1351  //================================================ Integration order
1352  this->determineIntegrationOrder ( maxDiag );
1353 
1354  //================================================ Report function completion
1355  ProSHADE_internal_messages::printProgressMessage ( this->verbose, 2, "Spherical harmonics environment prepared." );
1356 
1357  //================================================ Done
1358  return ;
1359 
1360 }

◆ determineBandwidth()

void ProSHADE_settings::determineBandwidth ( proshade_unsign  circumference)

This function determines the bandwidth for the spherical harmonics computation.

This function is here to automstically determine the bandwidth to which the spherical harmonics computations should be done. It accomplishes this by checking if value is already set, and if not (value is 0), then it sets it to half of the maximum circumference of the map, in indices as recommended by Kostelec and Rockmore (2007).

Parameters
[in]circumferenceThe maximum circumference of the map.

Definition at line 1180 of file ProSHADE.cpp.

1181 {
1182  //================================================ Check the current settings value is set to auto
1183  if ( this->maxBandwidth != 0 )
1184  {
1185  std::stringstream hlpSS;
1186  hlpSS << "The bandwidth was determined as: " << this->maxBandwidth;
1188  return ;
1189  }
1190 
1191  //================================================ Determine automatically
1193 
1194  //================================================ Report progress
1195  std::stringstream hlpSS;
1196  hlpSS << "The bandwidth was determined as: " << this->maxBandwidth;
1198 
1199  //================================================ Done
1200  return ;
1201 
1202 }

◆ determineBandwidthFromAngle()

void ProSHADE_settings::determineBandwidthFromAngle ( proshade_double  uncertainty)

This function determines the bandwidth for the spherical harmonics computation from the allowed rotation function angle uncertainty.

This function makes use of the fact that the rotation function dimensions will be 2 * bandwidth and that the dimensions will be covering full 360 degrees rotation space. Therefore, by saying what is the maximum allowed angle uncertainty, the minimum required bandwidth value can be determined.

Parameters
[in]uncertaintyThe maximum allowed uncertainty on the rotation function.

Definition at line 1212 of file ProSHADE.cpp.

1213 {
1214  //================================================ Determine bandwidth
1215  if ( static_cast<proshade_unsign> ( std::ceil ( ( 360.0 / uncertainty ) / 2 ) ) % 2 == 0 )
1216  {
1217  this->maxBandwidth = static_cast<proshade_unsign> ( std::ceil ( ( 360.0 / uncertainty ) / 2.0 ) );
1218  }
1219  else
1220  {
1221  this->maxBandwidth = static_cast<proshade_unsign> ( std::ceil ( ( 360.0 / uncertainty ) / 2.0 ) ) + 1;
1222  }
1223 
1224  //================================================ Report progress
1225  std::stringstream hlpSS;
1226  hlpSS << "The bandwidth was determined from uncertainty " << uncertainty << " degrees as: " << this->maxBandwidth;
1228 
1229  //================================================ Done
1230  return ;
1231 
1232 }

◆ determineIntegrationOrder()

void ProSHADE_settings::determineIntegrationOrder ( proshade_single  maxMapRange)

This function determines the integration order for the between spheres integration.

This function determines the order of the Gauss-Legendre integration which needs to be done between the spheres. To do this, it uses the pre-coputed values of maxium distance between integration points for each order and the maxium distance between spheres expressed as a fraction of the total.

Parameters
[in]maxMapRangeThe maximum diagonal distance of the map in Angstroms.

Definition at line 1275 of file ProSHADE.cpp.

1276 {
1277  //================================================ Check the current settings value is set to auto
1278  if ( this->integOrder != 0 )
1279  {
1280  std::stringstream hlpSS;
1281  hlpSS << "The integration order was determined as " << this->integOrder;
1283  return ;
1284  }
1285 
1286  //================================================ Determine automatically
1288 
1289  //================================================ Report progress
1290  std::stringstream hlpSS;
1291  hlpSS << "The integration order was determined as " << this->integOrder;
1293 
1294  //================================================ Done
1295  return ;
1296 
1297 }

◆ determineSphereDistances()

void ProSHADE_settings::determineSphereDistances ( proshade_single  maxMapRange)

This function determines the sphere distances for sphere mapping.

This function determines the distance between two consecutive spheres in the sphere mappin galgorithm. It checks if this values has not been already set and if not, it sets it as the sampling rate (distance between any two map points). It then checks that there will be at least 10 spheres and if not, it changes the sphere distance until at least 10 spheres are to be produced.

Parameters
[in]maxMapRangeThe maximum diagonal distance of the map in Angstroms.

Definition at line 1243 of file ProSHADE.cpp.

1244 {
1245  //================================================ Check the current settings value is set to auto
1246  if ( this->maxSphereDists != 0.0 )
1247  {
1248  std::stringstream hlpSS;
1249  hlpSS << "The sphere distances were determined as " << this->maxSphereDists << " Angstroms.";
1251  return ;
1252  }
1253 
1254  //================================================ Determine automatically
1256 
1257  //================================================ Report progress
1258  std::stringstream hlpSS;
1259  hlpSS << "The sphere distances were determined as " << this->maxSphereDists << " Angstroms.";
1261 
1262  //================================================ Done
1263  return ;
1264 
1265 }

◆ getCommandLineParams()

void ProSHADE_settings::getCommandLineParams ( int  argc,
char **  argv 
)

This function parses the command line arguments into the settings object.

Parameters
[in]argcThe count of the command line arguments (as passed to main function by the system).
[in]argvThe string containing the command line arguments (as passed to main function by the system).

Definition at line 1583 of file ProSHADE.cpp.

1584 {
1585  //================================================ If no command line arguments, print help
1586  if ( argc == 1 ) { ProSHADE_internal_messages::printHelp ( ); }
1587 
1588  //================================================ Long options struct
1589  const struct option_port longopts[] =
1590  {
1591  { "version", no_argument, NULL, 'v' },
1592  { "help", no_argument, NULL, 'h' },
1593  { "verbose", required_argument, NULL, '!' },
1594  { "distances", no_argument, NULL, 'D' },
1595  { "mapManip", no_argument, NULL, 'M' },
1596  { "symmetry", no_argument, NULL, 'S' },
1597  { "overlay", no_argument, NULL, 'O' },
1598  { "file", required_argument, NULL, 'f' },
1599  { "forceSpgP1", no_argument, NULL, 'u' },
1600  { "removeWaters", no_argument, NULL, 'w' },
1601  { "firstModel", no_argument, NULL, 'x' },
1602  { "resolution", required_argument, NULL, 'r' },
1603  { "bandwidth", required_argument, NULL, 'b' },
1604  { "sphereDists", required_argument, NULL, 's' },
1605  { "extraSpace", required_argument, NULL, 'e' },
1606  { "integOrder", required_argument, NULL, 'i' },
1607  { "taylorCap", required_argument, NULL, 't' },
1608  { "invertMap", no_argument, NULL, '@' },
1609  { "normalise", no_argument, NULL, '#' },
1610  { "mask", no_argument, NULL, '$' },
1611  { "saveMask", no_argument, NULL, '%' },
1612  { "maskFile", required_argument, NULL, '^' },
1613  { "maskBlurring", required_argument, NULL, '&' },
1614  { "maskThreshold", required_argument, NULL, '*' },
1615  { "mapReboxing", no_argument, NULL, 'R' },
1616  { "boundsSpace", required_argument, NULL, '(' },
1617  { "boundsThreshold", required_argument, NULL, ')' },
1618  { "sameBoundaries", no_argument, NULL, '-' },
1619  { "reBoxedFilename", required_argument, NULL, 'g' },
1620  { "pdbTempFact", required_argument, NULL, 'd' },
1621  { "center", no_argument, NULL, 'c' },
1622  { "changeMapResol", no_argument, NULL, 'j' },
1623  { "changeMapTriLin", no_argument, NULL, 'a' },
1624  { "noPhase", no_argument, NULL, 'p' },
1625  { "progressive", no_argument, NULL, 'k' },
1626  { "noEnL", no_argument, NULL, 'l' },
1627  { "noTrS", no_argument, NULL, 'm' },
1628  { "noFRF", no_argument, NULL, 'n' },
1629  { "EnLWeight", required_argument, NULL, '_' },
1630  { "peakNeigh", required_argument, NULL, '=' },
1631  { "peakThres", required_argument, NULL, '+' },
1632  { "missAxThres", required_argument, NULL, '[' },
1633  { "sameAxComp", required_argument, NULL, ']' },
1634  { "axisComBeh", no_argument, NULL, 'q' },
1635  { "bicubSearch", no_argument, NULL, 'A' },
1636  { "maxSymPrime", required_argument, NULL, 'B' },
1637  { "minPeakHeight", required_argument, NULL, 'o' },
1638  { "reqSym", required_argument, NULL, '{' },
1639  { "overlayFile", required_argument, NULL, '}' },
1640  { "overlayJSONFile", required_argument, NULL, 'y' },
1641  { "angUncertain", required_argument, NULL, ';' },
1642  { "usePeaksInRotFun",no_argument, NULL, 'z' },
1643  { NULL, 0, NULL, 0 }
1644  };
1645 
1646  //================================================ Short options string
1647  const char* const shortopts = "AaB:b:cd:De:f:g:hi:jklmMno:Opqr:Rs:St:uvwxy:z!:@#$%^:&:*:(:):-_:=:+:[:]:{:}:;:";
1648 
1649  //================================================ Parsing the options
1650  while ( true )
1651  {
1652  //============================================ Read the next option
1653  int opt = getopt_long_port ( argc, argv, shortopts, longopts, NULL );
1654 
1655  //============================================ Done parsing
1656  if ( opt == -1 )
1657  {
1658  break;
1659  }
1660 
1661  //============================================ For each option, set the internal values appropriately
1662  switch ( opt )
1663  {
1664  //======================================= Print version info
1665  case 'v':
1666  {
1668  exit ( EXIT_SUCCESS );
1669  }
1670 
1671  //======================================= User needs help
1672  case 'h':
1673  {
1675  exit ( EXIT_SUCCESS );
1676  }
1677 
1678  //======================================= Save the argument as the verbosity value, or if no value given, just set to 3
1679  case '!':
1680  {
1681  this->setVerbosity ( static_cast<proshade_single> ( atoi ( optarg ) ) );
1682  continue;
1683  }
1684 
1685  //======================================= Set task to distances
1686  case 'D':
1687  {
1688  this->task = Distances;
1689  continue;
1690  }
1691 
1692  //======================================= Set task to map manipulation
1693  case 'M':
1694  {
1695  this->task = MapManip;
1696  continue;
1697  }
1698 
1699  //======================================= Set task to symmetry detection
1700  case 'S':
1701  {
1702  this->task = Symmetry;
1703 
1704  //=================================== Force default unless changed already by the user
1705  if ( this->requestedResolution == -1 ) { this->requestedResolution = 6.0; }
1706  if ( this->pdbBFactorNewVal == -1 ) { this->pdbBFactorNewVal = 80.0; }
1707  this->changeMapResolution = !this->changeMapResolution; // Switch value. This can be over-ridden by the user by using -j
1708  this->moveToCOM = !this->moveToCOM; // Switch value. This can be over-ridden by the user by using -c.
1709 
1710  continue;
1711  }
1712 
1713  //======================================= Set task to map overlay
1714  case 'O':
1715  {
1716  this->task = OverlayMap;
1717  continue;
1718  }
1719 
1720  //======================================= Save the argument as a file to read in
1721  case 'f':
1722  {
1723  this->addStructure ( std::string ( optarg ) );
1724  continue;
1725  }
1726 
1727  //======================================= Force the input PDB files to have P1 spacegroup
1728  case 'u':
1729  {
1730  this->forceP1 = !this->forceP1;
1731  continue;
1732  }
1733 
1734  //======================================= Remove waters from PDB input files?
1735  case 'w':
1736  {
1737  this->removeWaters = !this->removeWaters;
1738  continue;
1739  }
1740 
1741  //======================================= Use all models, or just the first one?
1742  case 'x':
1743  {
1744  this->firstModelOnly = !this->firstModelOnly;
1745  continue;
1746  }
1747 
1748  //======================================= Save the argument as the resolution value
1749  case 'r':
1750  {
1751  this->setResolution ( static_cast<proshade_single> ( atof ( optarg ) ) );
1752  continue;
1753  }
1754 
1755  //======================================= Save the argument as the bandwidth value
1756  case 'b':
1757  {
1758  this->setBandwidth ( static_cast<proshade_unsign> ( atoi ( optarg ) ) );
1759  continue;
1760  }
1761 
1762  //======================================= Save the argument as the extra space value
1763  case 'e':
1764  {
1765  this->setExtraSpace ( static_cast<proshade_single> ( atof ( optarg ) ) );
1766  continue;
1767  }
1768 
1769  //======================================= Save the argument as the intaggration order value
1770  case 'i':
1771  {
1772  this->setIntegrationOrder ( static_cast<proshade_unsign> ( atof ( optarg ) ) );
1773  continue;
1774  }
1775 
1776  //======================================= Save the argument as the sphere distance value
1777  case 's':
1778  {
1779  this->setSphereDistances ( static_cast<proshade_single> ( atof ( optarg ) ) );
1780  continue;
1781  }
1782 
1783  //======================================= Save the argument as the taylor series cap value
1784  case 't':
1785  {
1786  this->setTaylorSeriesCap ( static_cast<proshade_unsign> ( atof ( optarg ) ) );
1787  continue;
1788  }
1789 
1790  //======================================= Set map inversion to true
1791  case '@':
1792  {
1793  this->setMapInversion ( true );
1794  continue;
1795  }
1796 
1797  //======================================= Set map normalisation to true
1798  case '#':
1799  {
1800  this->setNormalisation ( true );
1801  continue;
1802  }
1803 
1804  //======================================= Set map masking to true
1805  case '$':
1806  {
1807  this->setMasking ( true );
1808  continue;
1809  }
1810 
1811  //======================================= Set map masking to true and mask map saving to true as well
1812  case '%':
1813  {
1814  this->setMasking ( true );
1815  this->setMaskSaving ( true );
1816  continue;
1817  }
1818 
1819  //======================================= Save the argument as the mask filename value
1820  case '^':
1821  {
1822  this->setMaskFilename ( static_cast<std::string> ( optarg ) );
1823  continue;
1824  }
1825 
1826  //======================================= Save the argument as the mask blurring factor value
1827  case '&':
1828  {
1829  this->setMaskBlurFactor ( static_cast<proshade_single> ( atof ( optarg ) ) );
1830  continue;
1831  }
1832 
1833  //======================================= Save the argument as the mask threshold (IQR) value
1834  case '*':
1835  {
1836  this->setMaskIQR ( static_cast<proshade_single> ( atof ( optarg ) ) );
1837  continue;
1838  }
1839 
1840  //======================================= Set map reboxing to true
1841  case 'R':
1842  {
1843  this->setMasking ( true );
1844  this->setMapReboxing ( true );
1845  continue;
1846  }
1847 
1848  //======================================= Save the argument as the bounds extra space value
1849  case '(':
1850  {
1851  this->setBoundsSpace ( static_cast<proshade_single> ( atof ( optarg ) ) );
1852  continue;
1853  }
1854 
1855  //======================================= Save the argument as the bounds similarity threshold value
1856  case ')':
1857  {
1858  this->setBoundsThreshold ( static_cast<proshade_signed> ( atoi ( optarg ) ) );
1859  continue;
1860  }
1861 
1862  //======================================= Set same boundaries to true
1863  case '-':
1864  {
1865  this->setSameBoundaries ( true );
1866  continue;
1867  }
1868 
1869  //======================================= Save the argument as the re-boxed structure filename value
1870  case 'g':
1871  {
1872  this->setOutputFilename ( static_cast<std::string> ( optarg ) );
1873  continue;
1874  }
1875 
1876  //======================================= Save the argument as the PDB B-factor new constant value
1877  case 'd':
1878  {
1879  this->setPDBBFactor ( static_cast<proshade_single> ( atof ( optarg ) ) );
1880  continue;
1881  }
1882 
1883  //======================================= Set map centering to true
1884  case 'c':
1885  {
1886  this->moveToCOM = !this->moveToCOM;
1887  continue;
1888  }
1889 
1890  //======================================= Set map resolution change using Fourier transforms to true
1891  case 'j':
1892  {
1894  continue;
1895  }
1896 
1897  //======================================= Set map resolution change using real-space tri-linear interpolation to true
1898  case 'a':
1899  {
1900  this->setMapResolutionChangeTriLinear ( true );
1901  continue;
1902  }
1903 
1904  //======================================= Set map phase removal to true
1905  case 'p':
1906  {
1907  this->setPhaseUsage ( false );
1908  continue;
1909  }
1910 
1911  //======================================= Set progressive shell mapping to true
1912  case 'k':
1913  {
1914  this->setProgressiveSphereMapping ( true );
1915  continue;
1916  }
1917 
1918  //======================================= Set energy level descriptor computation to false
1919  case 'l':
1920  {
1921  this->setEnergyLevelsComputation ( false );
1922  continue;
1923  }
1924 
1925  //======================================= Set trace sigma descriptor computation to false
1926  case 'm':
1927  {
1928  this->setTraceSigmaComputation ( false );
1929  continue;
1930  }
1931 
1932  //======================================= Set full rotation function descriptor computation to false
1933  case 'n':
1934  {
1935  this->setRotationFunctionComputation ( false );
1936  continue;
1937  }
1938 
1939  //======================================= Save the argument as the energy levels descriptor weight value
1940  case '_':
1941  {
1942  this->setEnLevShellWeight ( static_cast<proshade_double> ( atof ( optarg ) ) );
1943  continue;
1944  }
1945 
1946  //======================================= Save the argument as the peak neighbours minimum value
1947  case '=':
1948  {
1949  this->setPeakNeighboursNumber ( static_cast<proshade_unsign> ( atoi ( optarg ) ) );
1950  continue;
1951  }
1952 
1953  //======================================= Save the argument as the peak IQR from median naive small peaks removal value
1954  case '+':
1955  {
1956  this->setPeakNaiveNoIQR ( static_cast<proshade_double> ( atof ( optarg ) ) );
1957  continue;
1958  }
1959 
1960  //======================================= Save the argument as the missing axis threshold value
1961  case '[':
1962  {
1963  this->setMissingPeakThreshold ( static_cast<proshade_double> ( atof ( optarg ) ) );
1964  continue;
1965  }
1966 
1967  //======================================= Save the argument as the missing axis threshold value
1968  case ']':
1969  {
1970  setAxisComparisonThreshold ( static_cast<proshade_double> ( atof ( optarg ) ) );
1971  continue;
1972  }
1973 
1974  //======================================= Save the argument as the missing axis threshold value
1975  case 'q':
1976  {
1977  setAxisComparisonThresholdBehaviour ( !this->axisErrToleranceDefault );
1978  continue;
1979  }
1980 
1981  //======================================= Save the argument as the bicubic interpolation search requirement value
1982  case 'A':
1983  {
1985  continue;
1986  }
1987 
1988  //======================================= Save the argument as the bicubic interpolation search requirement value
1989  case 'B':
1990  {
1991  setMaxSymmetryFold ( static_cast<proshade_unsign> ( atoi ( optarg ) ) );
1992  continue;
1993  }
1994 
1995  //======================================= Minimum peak height for axis
1996  case 'o':
1997  {
1998  this->minSymPeak = static_cast<proshade_double> ( atof ( optarg ) );
1999  continue;
2000  }
2001 
2002  //======================================= Save the argument as the requested symmetry and potentially fold value
2003  case '{':
2004  {
2005  std::string input = static_cast<std::string> ( optarg );
2006 
2007  if ( input.at(0) == 'C' )
2008  {
2009  this->setRequestedSymmetry ( "C" );
2010 
2011  std::string numHlp ( input.begin()+1, input.end() );
2012  if ( numHlp.length() > 0 ) { this->setRequestedFold ( atoi ( numHlp.c_str() ) ); }
2013  else { std::cerr << "!!! ProSHADE ERROR !!! The input argument requests search for Cyclic/Dihedral symmetry, but does not specify the requested fold." << std::endl; exit ( EXIT_FAILURE ); }
2014  }
2015  else
2016  {
2017  if ( input.at(0) == 'D' )
2018  {
2019  this->setRequestedSymmetry ( "D" );
2020 
2021  std::string numHlp ( input.begin()+1, input.end() );
2022  if ( numHlp.length() > 0 ) { this->setRequestedFold ( atoi ( numHlp.c_str() ) ); }
2023  else { std::cerr << "!!! ProSHADE ERROR !!! The input argument requests search for Cyclic/Dihedral symmetry, but does not specify the requested fold." << std::endl; exit ( EXIT_FAILURE ); }
2024  }
2025  else
2026  {
2027  if ( input.at(0) == 'T' )
2028  {
2029  this->setRequestedSymmetry ( "T" );
2030  }
2031  else
2032  {
2033  if ( input.at(0) == 'O' )
2034  {
2035  this->setRequestedSymmetry ( "O" );
2036  }
2037  else
2038  {
2039  if ( input.at(0) == 'I' )
2040  {
2041  this->setRequestedSymmetry ( "I" );
2042  }
2043  else
2044  {
2045  std::cerr << "!!! ProSHADE ERROR !!! Failed to parse the requested symmetry type. Allowed types are C, D, T, O and I, with C and D requiring to be followed by a number specifying the fold." << std::endl; exit ( EXIT_FAILURE );
2046  }
2047  }
2048  }
2049  }
2050  }
2051 
2052  continue;
2053  }
2054 
2055  //======================================= Save the argument as filename to save the overlay moved structure to value
2056  case '}':
2057  {
2058  this->setOverlaySaveFile ( static_cast<std::string> ( optarg ) );
2059  continue;
2060  }
2061 
2062  //======================================= Save the argument as filename to save the overlay operations to value
2063  case 'y':
2064  {
2065  this->setOverlayJsonFile ( static_cast<std::string> ( optarg ) );
2066  continue;
2067  }
2068 
2069  //======================================= Save the argument as angular uncertainty for bandwidth determination
2070  case ';':
2071  {
2072  this->rotationUncertainty = static_cast<proshade_double> ( atof ( optarg ) );
2073  continue;
2074  }
2075 
2076  //======================================= Save the argument as angular uncertainty for bandwidth determination
2077  case 'z':
2078  {
2079  this->setSymmetryRotFunPeaks ( false );
2080  continue;
2081  }
2082 
2083  //======================================= Unknown option
2084  case '?':
2085  {
2086  //=================================== Save the argument as angular uncertainty for bandwidth determination
2087  if ( optopt )
2088  {
2089  std::cout << "!!! ProSHADE ERROR !!! Unrecognised short option -" << static_cast<char> ( optopt ) << " . Please use -h for help on the command line options." << std::endl;
2090  }
2091  else
2092  {
2093  std::cout << "!!! ProSHADE ERROR !!! Unrecognised long option " << argv[static_cast<int> (optind)-1] << " . Please use -h for help on the command line options." << std::endl;
2094  }
2095 
2096  //=================================== This case is handled by getopt_long, nothing more needed.
2097  exit ( EXIT_SUCCESS );
2098  }
2099 
2100  //======================================= Fallback option
2101  default:
2102  {
2104  exit ( EXIT_SUCCESS );
2105  }
2106  }
2107  }
2108 
2109  //================================================ Done
2110  return ;
2111 
2112 }

◆ printSettings()

void ProSHADE_settings::printSettings ( void  )

This function prints the current values in the settings object.

Warning
This is a debugging function of no real utility to the user.

Definition at line 2118 of file ProSHADE.cpp.

2119 {
2120  //================================================ Print the currest values in the settings object
2121  std::stringstream strstr;
2122  strstr.str(std::string());
2123  if ( this->task == NA ) { strstr << "NA"; }
2124  if ( this->task == Distances ) { strstr << "DISTANCES COMPUTATION"; }
2125  if ( this->task == MapManip ) { strstr << "MAP MANIPULATION"; }
2126  if ( this->task == Symmetry ) { strstr << "SYMMETRY DETECTION"; }
2127  if ( this->task == OverlayMap ) { strstr << "MAP OVERLAY"; }
2128  printf ( "Task to perform : %37s\n", strstr.str().c_str() );
2129 
2130  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( this->inputFiles.size() ); iter++ )
2131  {
2132  strstr.str(std::string());
2133  strstr << this->inputFiles.at(iter);
2134  printf ( "File(s) to process : %37s\n", strstr.str().c_str() );
2135  }
2136 
2137  strstr.str(std::string());
2138  strstr << this->verbose;
2139  printf ( "Verbosity : %37s\n", strstr.str().c_str() );
2140 
2141  strstr.str(std::string());
2142  strstr << this->requestedResolution;
2143  printf ( "Resolution (comp) : %37s\n", strstr.str().c_str() );
2144 
2145  strstr.str(std::string());
2146  strstr << this->maxBandwidth;
2147  printf ( "Bandwidth : %37s\n", strstr.str().c_str() );
2148 
2149  strstr.str(std::string());
2150  strstr << this->maxSphereDists;
2151  printf ( "Sphere distances : %37s\n", strstr.str().c_str() );
2152 
2153  strstr.str(std::string());
2154  strstr << this->addExtraSpace;
2155  printf ( "Extra space : %37s\n", strstr.str().c_str() );
2156 
2157  strstr.str(std::string());
2158  if ( this->forceP1 ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2159  printf ( "Force P1 spacegroup : %37s\n", strstr.str().c_str() );
2160 
2161  strstr.str(std::string());
2162  if ( this->removeWaters ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2163  printf ( "Waters removed : %37s\n", strstr.str().c_str() );
2164 
2165  strstr.str(std::string());
2166  if ( this->firstModelOnly ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2167  printf ( "Only 1st model : %37s\n", strstr.str().c_str() );
2168 
2169  strstr.str(std::string());
2170  strstr << this->integOrder;
2171  printf ( "Integration order : %37s\n", strstr.str().c_str() );
2172 
2173  strstr.str(std::string());
2174  strstr << this->taylorSeriesCap;
2175  printf ( "Taylor series cap : %37s\n", strstr.str().c_str() );
2176 
2177  strstr.str(std::string());
2178  strstr << this->pdbBFactorNewVal;
2179  printf ( "PDB B-factor const : %37s\n", strstr.str().c_str() );
2180 
2181  strstr.str(std::string());
2182  if ( this->reBoxMap ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2183  printf ( "Map re-boxing : %37s\n", strstr.str().c_str() );
2184 
2185  strstr.str(std::string());
2186  if ( this->invertMap ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2187  printf ( "Map inversion : %37s\n", strstr.str().c_str() );
2188 
2189  strstr.str(std::string());
2190  if ( this->normaliseMap ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2191  printf ( "Map normalisation : %37s\n", strstr.str().c_str() );
2192 
2193  strstr.str(std::string());
2194  if ( this->maskMap ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2195  printf ( "Map masking : %37s\n", strstr.str().c_str() );
2196 
2197  strstr.str(std::string());
2198  if ( this->saveMask ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2199  printf ( "Saving mask : %37s\n", strstr.str().c_str() );
2200 
2201  strstr.str(std::string());
2202  strstr << this->maskFileName;
2203  printf ( "Map mask filename : %37s\n", strstr.str().c_str() );
2204 
2205  strstr.str(std::string());
2206  strstr << this->blurFactor;
2207  printf ( "Map blurring : %37s\n", strstr.str().c_str() );
2208 
2209  strstr.str(std::string());
2210  strstr << this->maskingThresholdIQRs;
2211  printf ( "Masking threshold : %37s\n", strstr.str().c_str() );
2212 
2213  strstr.str(std::string());
2214  strstr << this->boundsExtraSpace;
2215  printf ( "Bounds extra space : %37s\n", strstr.str().c_str() );
2216 
2217  strstr.str(std::string());
2218  strstr << this->boundsSimilarityThreshold;
2219  printf ( "Bounds similarity : %37s\n", strstr.str().c_str() );
2220 
2221  strstr.str(std::string());
2222  if ( this->useSameBounds ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2223  printf ( "Same boundaries : %37s\n", strstr.str().c_str() );
2224 
2225  strstr.str(std::string());
2226  strstr << this->outName;
2227  printf ( "Re-boxed filename : %37s\n", strstr.str().c_str() );
2228 
2229  strstr.str(std::string());
2230  if ( this->moveToCOM ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2231  printf ( "Map COM centering : %37s\n", strstr.str().c_str() );
2232 
2233  strstr.str(std::string());
2234  if ( this->changeMapResolution ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2235  printf ( "Change map resol : %37s\n", strstr.str().c_str() );
2236 
2237  strstr.str(std::string());
2238  if ( this->changeMapResolutionTriLinear ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2239  printf ( "Change map tri-lin : %37s\n", strstr.str().c_str() );
2240 
2241  strstr.str(std::string());
2242  if ( this->usePhase ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2243  printf ( "Use phase info : %37s\n", strstr.str().c_str() );
2244 
2245  strstr.str(std::string());
2246  if ( this->progressiveSphereMapping ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2247  printf ( "Progressive spheres : %37s\n", strstr.str().c_str() );
2248 
2249  strstr.str(std::string());
2250  if ( this->computeEnergyLevelsDesc ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2251  printf ( "Energy lvl desc : %37s\n", strstr.str().c_str() );
2252 
2253  strstr.str(std::string());
2254  strstr << this->enLevMatrixPowerWeight;
2255  printf ( "Energy lvl weight : %37s\n", strstr.str().c_str() );
2256 
2257  strstr.str(std::string());
2258  if ( this->computeTraceSigmaDesc ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2259  printf ( "Tr sigma desc : %37s\n", strstr.str().c_str() );
2260 
2261  strstr.str(std::string());
2262  if ( this->computeRotationFuncDesc ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2263  printf ( "Full RF desc : %37s\n", strstr.str().c_str() );
2264 
2265  strstr.str(std::string());
2266  strstr << this->peakNeighbours;
2267  printf ( "Neightbours to peak : %37s\n", strstr.str().c_str() );
2268 
2269  strstr.str(std::string());
2270  strstr << this->noIQRsFromMedianNaivePeak;
2271  printf ( "Peak IQR threshold : %37s\n", strstr.str().c_str() );
2272 
2273  strstr.str(std::string());
2274  strstr << this->symMissPeakThres;
2275  printf ( "Missing ax. thres : %37s\n", strstr.str().c_str() );
2276 
2277  strstr.str(std::string());
2278  strstr << this->minSymPeak;
2279  printf ( "Min. sym. peak size : %37s\n", strstr.str().c_str() );
2280 
2281  strstr.str(std::string());
2282  strstr << this->axisErrTolerance;
2283  printf ( "Same ax. threshold : %37s\n", strstr.str().c_str() );
2284 
2285  strstr.str(std::string());
2286  if ( this->axisErrToleranceDefault ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2287  printf ( "Same ax. thre. decr.: %37s\n", strstr.str().c_str() );
2288 
2289  strstr.str(std::string());
2290  strstr << this->requestedSymmetryType << "-" << this->requestedSymmetryFold;
2291  printf ( "Requested symm. : %37s\n", strstr.str().c_str() );
2292 
2293  strstr.str(std::string());
2294  strstr << this->overlayStructureName;
2295  printf ( "Overlay file : %37s\n", strstr.str().c_str() );
2296 
2297  strstr.str(std::string());
2298  strstr << this->rotTrsJSONFile;
2299  printf ( "JSON overlay file : %37s\n", strstr.str().c_str() );
2300 
2301  //================================================ Done
2302  return ;
2303 
2304 }

◆ setAxisComparisonThreshold()

void ProSHADE_settings::setAxisComparisonThreshold ( proshade_double  axThres)

Sets the threshold for matching symmetry axes.

When comparing symmetry axes, there needs to be a threshold allowing for some small error comming from the numberical inaccuracies. This is where you set this threshold.

Parameters
[in]axThresThe requested value for the axes comparison threshold.

Definition at line 955 of file ProSHADE.cpp.

956 {
957  //================================================ Set the value
958  this->axisErrTolerance = axThres;
959 
960  //================================================ Done
961  return ;
962 
963 }

◆ setAxisComparisonThresholdBehaviour()

void ProSHADE_settings::setAxisComparisonThresholdBehaviour ( bool  behav)

Sets the automatic symmetry axis tolerance decreasing.

When comparing symmetry axes, there needs to be a threshold allowing for some small error comming from the numberical inaccuracies. It turns out that this threshold should take into account the ratio to the next symmetry angles, otherwise it would strongly prefer larger symmetries. This variable decides whether the threshold should be decreased based on the fold of sought åsymmetry or not.

Parameters
[in]behavThe requested value for the axes comparison threshold decreasing.

Definition at line 974 of file ProSHADE.cpp.

975 {
976  //================================================ Set the value
977  this->axisErrToleranceDefault = behav;
978 
979  //================================================ Done
980  return ;
981 
982 }

◆ setBandwidth()

void ProSHADE_settings::setBandwidth ( proshade_unsign  band)

Sets the requested spherical harmonics bandwidth in the appropriate variable.

This function sets the bandwidth limit for the spherical harmonics computations in the appropriate variable.

Parameters
[in]bandThe requested value for spherical harmonics bandwidth (0 = AUTOMATIC DETERMINATION).

Definition at line 732 of file ProSHADE.cpp.

733 {
734  //================================================ Set the value
735  this->maxBandwidth = band;
736 
737  //================================================ Done
738  return ;
739 
740 }

◆ setBicubicInterpolationSearch()

void ProSHADE_settings::setBicubicInterpolationSearch ( bool  bicubPeaks)

Sets the bicubic interpolation on peaks.

Parameters
[in]bicubPeaksShould bicubic interpolation be done to search for improved axis in between peak index values (DEFAULT - TRUE)?

Definition at line 1148 of file ProSHADE.cpp.

1149 {
1150  //================================================ Set the value
1151  this->useBiCubicInterpolationOnPeaks = bicubPeaks;
1152 
1153  //================================================ Done
1154  return ;
1155 
1156 }

◆ setBoundsSpace()

void ProSHADE_settings::setBoundsSpace ( proshade_single  boundsExSp)

Sets the requested number of angstroms for extra space in re-boxing in the appropriate variable.

This function sets the number of angstroms to be added both before and after the absolute bounds for re-boxing to the correct variable.

Parameters
[in]boundsExSpThe requested value for the extra re-boxing space in angstroms.

Definition at line 586 of file ProSHADE.cpp.

587 {
588  //================================================ Set the value
589  this->boundsExtraSpace = boundsExSp;
590 
591  //================================================ Done
592  return ;
593 
594 }

◆ setBoundsThreshold()

void ProSHADE_settings::setBoundsThreshold ( proshade_signed  boundsThres)

Sets the threshold for number of indices difference acceptable to make index sizes same in the appropriate variable.

This function sets the number of indices by which two dimensions can differ for them to be still made the same size.

Parameters
[in]boundsThresThe requested value for the bouds difference threhshold.

Definition at line 602 of file ProSHADE.cpp.

603 {
604  //================================================ Set the value
605  this->boundsSimilarityThreshold = boundsThres;
606 
607  //================================================ Done
608  return ;
609 
610 }

◆ setCorrelationMasking()

void ProSHADE_settings::setCorrelationMasking ( bool  corMask)

Sets the requested map masking type in the appropriate variable.

This function sets the map masking type. If false, the standard map blurring masking will be used, while if true, the new "fake" half-map correlation mask will be used.

Parameters
[in]corMaskThe requested value for the map masking type.

Definition at line 487 of file ProSHADE.cpp.

488 {
489  //================================================ Set the value
490  this->useCorrelationMasking = corMask;
491 
492  //================================================ Done
493  return ;
494 
495 }

◆ setDetectedSymmetry()

void ProSHADE_settings::setDetectedSymmetry ( proshade_double *  sym)

Sets the final detected symmetry axes information.

This function copies (deep copy) the detected and recommended (or requested) symmetry axis information into the settings object variable for further processing. For multiple axes, call this function multiple times - the addition is cumulative.

Parameters
[in]symA pointer to single symmetry axis constituting the detected symmetry.

Definition at line 1077 of file ProSHADE.cpp.

1078 {
1079  //================================================ Allocate memory
1080  proshade_double* hlpAxis = new proshade_double [6];
1081  ProSHADE_internal_misc::checkMemoryAllocation ( hlpAxis, __FILE__, __LINE__, __func__ );
1082 
1083  //================================================ Copy (deep) data
1084  hlpAxis[0] = sym[0];
1085  hlpAxis[1] = sym[1];
1086  hlpAxis[2] = sym[2];
1087  hlpAxis[3] = sym[3];
1088  hlpAxis[4] = sym[4];
1089  hlpAxis[5] = sym[5];
1090 
1091  //================================================ Save
1093 
1094  //================================================ Release memory
1095  delete[] hlpAxis;
1096 
1097  //================================================ Done
1098  return ;
1099 
1100 }

◆ setEnergyLevelsComputation()

void ProSHADE_settings::setEnergyLevelsComputation ( bool  enLevDesc)

Sets whether the energy level distance descriptor should be computed.

This function sets the boolean variable deciding whether the RRP matrices and the energy levels descriptor should be computed or not.

Parameters
[in]enLevDescThe requested value for the energy levels descriptor computation switch.

Definition at line 798 of file ProSHADE.cpp.

799 {
800  //======================================== Set the value
801  this->computeEnergyLevelsDesc = enLevDesc;
802 
803  //======================================== Done
804  return ;
805 
806 }

◆ setEnLevShellWeight()

void ProSHADE_settings::setEnLevShellWeight ( proshade_double  mPower)

Sets the weight of shell position for the energy levels computation.

During the computation of the energy levels descriptor, Pearson's correlation coefficient is computed between different shells with the same band. The shell index can by expanded to its mPower exponential to give higher shells more weight, or vice versa. To do this, set the mPower value as you see fit.

Parameters
[in]mPowerThe requested value for the shell position exponential.

Definition at line 903 of file ProSHADE.cpp.

904 {
905  //================================================ Set the value
906  this->enLevMatrixPowerWeight = mPower;
907 
908  //================================================ Done
909  return ;
910 
911 }

◆ setExtraSpace()

void ProSHADE_settings::setExtraSpace ( proshade_single  exSpace)

Sets the requested map extra space value in the appropriate variable.

This function sets the amount of extra space to be added to internal maps in the appropriate variable.

Parameters
[in]exSpaceThe requested amount of extra space.

Definition at line 700 of file ProSHADE.cpp.

701 {
702  //================================================ Set the value
703  this->addExtraSpace = exSpace;
704 
705  //================================================ Done
706  return ;
707 
708 }

◆ setGroupingSmoothingFactor()

void ProSHADE_settings::setGroupingSmoothingFactor ( proshade_double  smFact)

Sets the grouping smoothing factor into the proper variable.

When detecting symmetry, it is worth grouping the possible rotations by their self-rotation function peak heights. In this process, the distribution of peak heights needs to be smoothen over and this factor decides how smooth it should be. Small value leads to all peaks being in the same group, while large number means each peak will be in its own group.

Parameters
[in]smFactThe requested value for the grouping smoothing factor.

Definition at line 921 of file ProSHADE.cpp.

922 {
923  //================================================ Set the value
924  this->smoothingFactor = smFact;
925 
926  //================================================ Done
927  return ;
928 
929 }

◆ setIntegrationOrder()

void ProSHADE_settings::setIntegrationOrder ( proshade_unsign  intOrd)

Sets the requested order for the Gauss-Legendre integration in the appropriate variable.

This function sets the Gauss-Legendre integration between the spheres order value in the appropriate variable.

Parameters
[in]intOrdThe requested value for the Gauss-Legendre integration order (0 = AUTOMATIC DETERMINATION).

Definition at line 764 of file ProSHADE.cpp.

765 {
766  //================================================ Set the value
767  this->integOrder = intOrd;
768 
769  //================================================ Done
770  return ;
771 
772 }

◆ setMapCentering()

void ProSHADE_settings::setMapCentering ( bool  com)

Sets the requested map centering decision value in the appropriate variable.

This function sets the map centering using COM between on and off.

Parameters
[in]comThe requested value for the map centering (on = true, off = false).

Definition at line 684 of file ProSHADE.cpp.

685 {
686  //================================================ Set the value
687  this->moveToCOM = com;
688 
689  //================================================ Done
690  return ;
691 
692 }

◆ setMapInversion()

void ProSHADE_settings::setMapInversion ( bool  mInv)

Sets the requested map inversion value in the appropriate variable.

This function sets the map inversion between on and off.

Parameters
[in]mInvShould the map be inverted? (on = true, off = false).

Definition at line 405 of file ProSHADE.cpp.

406 {
407  //================================================ Set the value
408  this->invertMap = mInv;
409 
410  //================================================ Done
411  return ;
412 
413 }

◆ setMapReboxing()

void ProSHADE_settings::setMapReboxing ( bool  reBx)

Sets whether re-boxing needs to be done in the appropriate variable.

This function sets the switch as to whether re-boxing needs to be done to the correct variable.

Parameters
[in]reBxThe requested value for the re-boxing switch variable.

Definition at line 569 of file ProSHADE.cpp.

570 {
571  //================================================ Set the value
572  this->reBoxMap = reBx;
573 
574  //================================================ Done
575  return ;
576 
577 }

◆ setMapResolutionChange()

void ProSHADE_settings::setMapResolutionChange ( bool  mrChange)

Sets the requested map resolution change decision in the appropriate variable.

This function sets the map resolution change between on and off.

Parameters
[in]mrChangeThe requested value for the map resolution change (on = true, off = false).

Definition at line 652 of file ProSHADE.cpp.

653 {
654  //================================================ Set the value
655  this->changeMapResolution = mrChange;
656 
657  //================================================ Done
658  return ;
659 
660 }

◆ setMapResolutionChangeTriLinear()

void ProSHADE_settings::setMapResolutionChangeTriLinear ( bool  mrChange)

Sets the requested map resolution change decision using tri-linear interpolation in the appropriate variable.

This function sets the tri-linear interpolation map resolution change between on and off.

Parameters
[in]mrChangeThe requested value for the map resolution change (on = true, off = false).

Definition at line 668 of file ProSHADE.cpp.

669 {
670  //================================================ Set the value
671  this->changeMapResolutionTriLinear = mrChange;
672 
673  //================================================ Done
674  return ;
675 
676 }

◆ setMaskBlurFactor()

void ProSHADE_settings::setMaskBlurFactor ( proshade_single  blurFac)

Sets the requested map blurring factor in the appropriate variable.

This function sets the blurring / sharpening factor for map masking in the appropriate variable.

Parameters
[in]blurFacThe requested value for the blurring factor.

Definition at line 437 of file ProSHADE.cpp.

438 {
439  //================================================ Set the value
440  this->blurFactor = blurFac;
441 
442  //================================================ Done
443  return ;
444 
445 }

◆ setMaskFilename()

void ProSHADE_settings::setMaskFilename ( std::string  mskFln)

Sets where the mask should be saved.

This function sets the the filename to which mask should be saved.

Parameters
[in]mskFlnThe filename where the mask should be saved.

Definition at line 553 of file ProSHADE.cpp.

554 {
555  //================================================ Set the value
556  this->maskFileName = mskFln;
557 
558  //================================================ Done
559  return ;
560 
561 }

◆ setMasking()

void ProSHADE_settings::setMasking ( bool  mask)

Sets the requested map masking decision value in the appropriate variable.

This function sets the map masking between on and off.

Parameters
[in]maskThe requested value for the map masking (on = true, off = false).

Definition at line 470 of file ProSHADE.cpp.

471 {
472  //================================================ Set the value
473  this->maskMap = mask;
474 
475  //================================================ Done
476  return ;
477 
478 }

◆ setMaskIQR()

void ProSHADE_settings::setMaskIQR ( proshade_single  noIQRs)

Sets the requested number of IQRs for masking threshold in the appropriate variable.

This function sets the number of interquartile ranges from the median to be used for map masking in the correct variable.

Parameters
[in]noIQRsThe requested value for the number of IQRs from the median to be used for masking threshold.

Definition at line 454 of file ProSHADE.cpp.

455 {
456  //================================================ Set the value
457  this->maskingThresholdIQRs = noIQRs;
458 
459  //================================================ Done
460  return ;
461 
462 }

◆ setMaskSaving()

void ProSHADE_settings::setMaskSaving ( bool  savMsk)

Sets whether the mask should be saved.

This function sets the switch variable to whether mask should be saved.

Parameters
[in]savMskIf true, mask will be saved, otherwise it will not be.

Definition at line 537 of file ProSHADE.cpp.

538 {
539  //================================================ Set the value
540  this->saveMask = savMsk;
541 
542  //================================================ Done
543  return ;
544 
545 }

◆ setMaxSymmetryFold()

void ProSHADE_settings::setMaxSymmetryFold ( proshade_unsign  maxFold)

Sets the maximum symmetry fold (well, the maximum prime symmetry fold).

Parameters
[in]maxFoldMaximum prime number fold that will be searched for. Still its multiples may also be found.

Definition at line 1162 of file ProSHADE.cpp.

1163 {
1164  //================================================ Set the value
1165  this->maxSymmetryFold = maxFold;
1166 
1167  //================================================ Done
1168  return ;
1169 
1170 }

◆ setMinimumMaskSize()

void ProSHADE_settings::setMinimumMaskSize ( proshade_single  minMS)

Sets the requested minimum mask size.

This function sets the kernel for the local correlation computation between the "fake half-map" and the original map.

Parameters
[in]minMSThe requested value for the minimum mask size in Angstrom.

Definition at line 521 of file ProSHADE.cpp.

522 {
523  //================================================ Set the value
524  this->correlationKernel = minMS;
525 
526  //================================================ Done
527  return ;
528 
529 }

◆ setMinimumPeakForAxis()

void ProSHADE_settings::setMinimumPeakForAxis ( proshade_double  minSP)

Sets the minimum peak height for symmetry axis to be considered.

When considering if a symmetry axis is "real" and should be acted upon, its average peak height will need to be higher than this value.

Parameters
[in]minSPThe requested value for the minimum peak height.

Definition at line 991 of file ProSHADE.cpp.

992 {
993  //================================================ Set the value
994  this->minSymPeak = minSP;
995 
996  //================================================ Done
997  return ;
998 
999 }

◆ setMissingPeakThreshold()

void ProSHADE_settings::setMissingPeakThreshold ( proshade_double  mpThres)

Sets the threshold for starting the missing peaks procedure.

When only mpThres percentage of peaks are missing during symmetry detection, the full missing peak detection procedure will be started. Otherwise, the symmetry will not be detected at all.

Parameters
[in]mpThresThe requested value for the missing peaks procedure starting threshold.

Definition at line 938 of file ProSHADE.cpp.

939 {
940  //================================================ Set the value
941  this->symMissPeakThres = mpThres;
942 
943  //================================================ Done
944  return ;
945 
946 }

◆ setNormalisation()

void ProSHADE_settings::setNormalisation ( bool  normalise)

Sets the requested map normalisation value in the appropriate variable.

This function sets the map normalisation between on and off.

Parameters
[in]normaliseThe requested value for the map normalisation (on = true, off = false).

Definition at line 389 of file ProSHADE.cpp.

390 {
391  //================================================ Set the value
392  this->normaliseMap = normalise;
393 
394  //================================================ Done
395  return ;
396 
397 }

◆ setOutputFilename()

void ProSHADE_settings::setOutputFilename ( std::string  oFileName)

Sets the requested output file name in the appropriate variable.

This function sets the filename to which the output structure(s) should be saved. This variable is used by multiple tasks and therefore cannot be more specifically described here.

Parameters
[in]oFileNameThe requested value for the output file name variable.

Definition at line 636 of file ProSHADE.cpp.

637 {
638  //================================================ Set the value
639  this->outName = oFileName;
640 
641  //================================================ Done
642  return ;
643 
644 }

◆ setOverlayJsonFile()

void ProSHADE_settings::setOverlayJsonFile ( std::string  filename)

Sets the filename to which the overlay operations are to be save into.

Parameters
[in]filenameThe filename to which the overlay operations are to be saved to.

Definition at line 1120 of file ProSHADE.cpp.

1121 {
1122  //================================================ Set the value
1123  this->rotTrsJSONFile = filename;
1124 
1125  //================================================ Done
1126  return ;
1127 
1128 }

◆ setOverlaySaveFile()

void ProSHADE_settings::setOverlaySaveFile ( std::string  filename)

Sets the filename to which the overlay structure is to be save into.

Parameters
[in]filenameThe filename to which the overlay structure is to be saved to.

Definition at line 1106 of file ProSHADE.cpp.

1107 {
1108  //================================================ Set the value
1109  this->overlayStructureName = filename;
1110 
1111  //================================================ Done
1112  return ;
1113 
1114 }

◆ setPDBBFactor()

void ProSHADE_settings::setPDBBFactor ( proshade_double  newBF)

Sets the requested B-factor value for PDB files in the appropriate variable.

This function sets the B-factor value for PDB files in the appropriate variable.

Parameters
[in]newBFThe requested value for the B-factor value for PDB files for smooth and processible maps.

Definition at line 373 of file ProSHADE.cpp.

374 {
375  //================================================ Set the value
376  this->pdbBFactorNewVal = newBF;
377 
378  //================================================ Done
379  return ;
380 
381 }

◆ setPeakNaiveNoIQR()

void ProSHADE_settings::setPeakNaiveNoIQR ( proshade_double  noIQRs)

Sets the number of IQRs from the median for threshold height a peak needs to be considered a peak.

This function sets the number of IQRs from the median that determine the threshold used to determine if a 'naive' peak is a peak, or just a random local maxim in the background. The set from which median and IQR is computed is the non-peak values.

Parameters
[in]noIQRsThe requested number of IQRs from the median.

Definition at line 867 of file ProSHADE.cpp.

868 {
869  //================================================ Set the value
870  this->noIQRsFromMedianNaivePeak = noIQRs;
871 
872  //================================================ Done
873  return ;
874 
875 }

◆ setPeakNeighboursNumber()

void ProSHADE_settings::setPeakNeighboursNumber ( proshade_unsign  pkS)

Sets the number of neighbour values that have to be smaller for an index to be considered a peak.

This function sets the number of neighbouring points (in all three dimensions and both positive and negative direction) that have to have lower value than the currently considered index in order for this index to be considered as a peak.

Parameters
[in]pkSThe requested value for the number of neighbours being lower for a peak.

Definition at line 849 of file ProSHADE.cpp.

850 {
851  //================================================ Set the value
852  this->peakNeighbours = pkS;
853 
854  //================================================ Done
855  return ;
856 
857 }

◆ setPhaseUsage()

void ProSHADE_settings::setPhaseUsage ( bool  phaseUsage)

Sets whether the phase information will be used.

This function sets the boolean variable deciding whether the phase information should be used. If not, Patterson maps will be used instead of density maps and the 3D data will be converted to them. Also, only even bands of the spherical harmonics decomposition will be computed as the odd bands must be 0.

Parameters
[in]phaseUsageThe requested value for the phase usage switch.

Definition at line 885 of file ProSHADE.cpp.

886 {
887  //================================================ Set the value
888  this->usePhase = phaseUsage;
889 
890  //================================================ Done
891  return ;
892 
893 }

◆ setProgressiveSphereMapping()

void ProSHADE_settings::setProgressiveSphereMapping ( bool  progSphMap)

Sets the requested sphere mapping value settings approach in the appropriate variable.

This function sets the progressive sphere mapping approach between on and off.

Parameters
[in]comThe requested value for the progressive sphere mapping (on = true, off = false).

Definition at line 716 of file ProSHADE.cpp.

717 {
718  //================================================ Set the value
719  this->progressiveSphereMapping = progSphMap;
720 
721  //================================================ Done
722  return ;
723 
724 }

◆ setRecommendedFold()

void ProSHADE_settings::setRecommendedFold ( proshade_unsign  val)

Sets the ProSHADE detected symmetry fold.

When symmetry detection is done, the resulting recommended symmetry fold (valid only for C and D symmetry types) will be saved in the settings object by this function.

Parameters
[in]valThe recommended symmetry fold for the structure.
Warning
This is an internal function and it should not be used by the user.

Definition at line 1028 of file ProSHADE.cpp.

1029 {
1030  //================================================ Set the value
1031  this->recommendedSymmetryFold = val;
1032 
1033  //================================================ Done
1034  return ;
1035 
1036 }

◆ setRecommendedSymmetry()

void ProSHADE_settings::setRecommendedSymmetry ( std::string  val)

Sets the ProSHADE detected symmetry type.

When symmetry detection is done, the resulting recommended symmetry type will be saved in the settings object by this function.

Parameters
[in]valThe recommended symmetry type for the structure.
Warning
This is an internal function and it should not be used by the user.

Definition at line 1009 of file ProSHADE.cpp.

1010 {
1011  //================================================ Set the value
1012  this->recommendedSymmetryType = val;
1013 
1014  //================================================ Done
1015  return ;
1016 
1017 }

◆ setRequestedFold()

void ProSHADE_settings::setRequestedFold ( proshade_unsign  val)

Sets the user requested symmetry fold.

When symmetry detection is started, this symmetry fold will be exclusively sought.

Parameters
[in]valThe requested symmetry fold for the structure.

Definition at line 1060 of file ProSHADE.cpp.

1061 {
1062  //================================================ Set the value
1063  this->requestedSymmetryFold = val;
1064 
1065  //================================================ Done
1066  return ;
1067 
1068 }

◆ setRequestedSymmetry()

void ProSHADE_settings::setRequestedSymmetry ( std::string  val)

Sets the user requested symmetry type.

When symmetry detection is started, this symmetry type will be exclusively sought.

Parameters
[in]valThe requested symmetry type for the structure.

Definition at line 1044 of file ProSHADE.cpp.

1045 {
1046  //================================================ Set the value
1047  this->requestedSymmetryType = val;
1048 
1049  //================================================ Done
1050  return ;
1051 
1052 }

◆ setResolution()

void ProSHADE_settings::setResolution ( proshade_single  resolution)

Sets the requested resolution in the appropriate variable.

This function sets the resolution in the appropriate variable.

Parameters
[in]resolutionThe requested value for the resolution to which the computations are to be done.

Definition at line 357 of file ProSHADE.cpp.

358 {
359  //================================================ Set the value
360  this->requestedResolution = resolution;
361 
362  //================================================ Done
363  return ;
364 
365 }

◆ setRotationFunctionComputation()

void ProSHADE_settings::setRotationFunctionComputation ( bool  rotfVal)

Sets whether the rotation function distance descriptor should be computed.

This function sets the boolean variable deciding whether the inverse SO(3) transform and the rotation function descriptor should be computed or not.

Parameters
[in]rotfValThe requested value for the rotation function descriptor computation switch.

Definition at line 832 of file ProSHADE.cpp.

833 {
834  //================================================ Set the value
835  this->computeRotationFuncDesc = rotfVal;
836 
837  //================================================ Done
838  return ;
839 
840 }

◆ setSameBoundaries()

void ProSHADE_settings::setSameBoundaries ( bool  sameB)

Sets whether same boundaries should be used in the appropriate variable.

This function sets the switch as to whether the same boundaries as for the first map should be forced upon the rest if the input maps.

Parameters
[in]sameBThe requested value for the same boundaries as first structure switch variable.

Definition at line 619 of file ProSHADE.cpp.

620 {
621  //================================================ Set the value
622  this->useSameBounds = sameB;
623 
624  //================================================ Done
625  return ;
626 
627 }

◆ setSphereDistances()

void ProSHADE_settings::setSphereDistances ( proshade_single  sphDist)

Sets the requested distance between spheres in the appropriate variable.

This function sets the distance between any two consecutive spheres in the sphere mapping of a map in the appropriate variable.

Parameters
[in]sphDistThe requested value for distance between spheres (0 = AUTOMATIC DETERMINATION).

Definition at line 748 of file ProSHADE.cpp.

749 {
750  //================================================ Set the value
751  this->maxSphereDists = sphDist;
752 
753  //================================================ Done
754  return ;
755 
756 }

◆ setSymmetryRotFunPeaks()

void ProSHADE_settings::setSymmetryRotFunPeaks ( bool  rotFunPeaks)

Sets the symmetry detection algorithm type.

Parameters
[in]rotFunPeaksShould the original peak detection in rotation function space be used (FALSE), or should the new angle-axis space search be used (DEFAULT - TRUE)?

Definition at line 1134 of file ProSHADE.cpp.

1135 {
1136  //================================================ Set the value
1137  this->usePeakSearchInRotationFunctionSpace = rotFunPeaks;
1138 
1139  //================================================ Done
1140  return ;
1141 
1142 }

◆ setTaylorSeriesCap()

void ProSHADE_settings::setTaylorSeriesCap ( proshade_unsign  tayCap)

Sets the requested Taylor series cap for the Gauss-Legendre integration in the appropriate variable.

This function sets the Taylor series maximum limit for the Gauss-Legendre integration between the spheres order value in the appropriate variable.

Parameters
[in]tayCapThe requested value for the Taylor series cap. (0 = AUTOMATIC DETERMINATION).

Definition at line 781 of file ProSHADE.cpp.

782 {
783  //================================================ Set the value
784  this->taylorSeriesCap = tayCap;
785 
786  //================================================ Done
787  return ;
788 
789 }

◆ setTraceSigmaComputation()

void ProSHADE_settings::setTraceSigmaComputation ( bool  trSigVal)

Sets whether the trace sigma distance descriptor should be computed.

This function sets the boolean variable deciding whether the E matrices and the trace sigma descriptor should be computed or not.

Parameters
[in]trSigValThe requested value for the trace sigma descriptor computation switch.

Definition at line 815 of file ProSHADE.cpp.

816 {
817  //================================================ Set the value
818  this->computeTraceSigmaDesc = trSigVal;
819 
820  //================================================ Done
821  return ;
822 
823 }

◆ setTypicalNoiseSize()

void ProSHADE_settings::setTypicalNoiseSize ( proshade_single  typNoi)

Sets the requested "fake" half-map kernel in the appropriate variable.

This function sets the kernel for creating the "fake" half-map. What is meant here is that a new map is created as the average of neighbours from the original map - this is useful in masking. This value then sets how many neighbours.

Parameters
[in]typNoiThe requested value for the typical noise size in Angstrom.

Definition at line 505 of file ProSHADE.cpp.

506 {
507  //================================================ Set the value
508  this->halfMapKernel = typNoi;
509 
510  //================================================ Done
511  return ;
512 
513 }

◆ setVerbosity()

void ProSHADE_settings::setVerbosity ( proshade_signed  verbosity)

Sets the requested verbosity in the appropriate variable.

This function sets the varbosity of the ProSHADE run in the appropriate variable.

Parameters
[in]verboseThe requested value for verbosity. -1 means no output, while 4 means loud output

Definition at line 421 of file ProSHADE.cpp.

422 {
423  //================================================ Set the value
424  this->verbose = verbosity;
425 
426  //================================================ Done
427  return ;
428 
429 }

The documentation for this class was generated from the following files:
ProSHADE_settings::noIQRsFromMedianNaivePeak
proshade_double noIQRsFromMedianNaivePeak
When doing peak searching, how many IQRs from the median the threshold for peak height should be (in ...
Definition: ProSHADE_settings.hpp:165
ProSHADE_settings::setOverlayJsonFile
void setOverlayJsonFile(std::string filename)
Sets the filename to which the overlay operations are to be save into.
Definition: ProSHADE.cpp:1120
ProSHADE_settings::maxBandwidth
proshade_unsign maxBandwidth
The bandwidth of spherical harmonics decomposition for the largest sphere.
Definition: ProSHADE_settings.hpp:109
ProSHADE_settings::integOrder
proshade_unsign integOrder
The order required for full Gauss-Legendre integration between the spheres.
Definition: ProSHADE_settings.hpp:119
ProSHADE_settings::recommendedSymmetryType
std::string recommendedSymmetryType
The symmetry type that ProSHADE finds the best fitting for the structure. Possible values are "" for ...
Definition: ProSHADE_settings.hpp:175
ProSHADE_settings::rotTrsJSONFile
std::string rotTrsJSONFile
The filename to which the rotation and translation operations are to be saved into.
Definition: ProSHADE_settings.hpp:185
ProSHADE_settings::setEnLevShellWeight
void setEnLevShellWeight(proshade_double mPower)
Sets the weight of shell position for the energy levels computation.
Definition: ProSHADE.cpp:903
ProSHADE_settings::computeTraceSigmaDesc
bool computeTraceSigmaDesc
If true, the trace sigma descriptor will be computed, otherwise all its computations will be omitted.
Definition: ProSHADE_settings.hpp:160
ProSHADE_settings::setTraceSigmaComputation
void setTraceSigmaComputation(bool trSigVal)
Sets whether the trace sigma distance descriptor should be computed.
Definition: ProSHADE.cpp:815
ProSHADE_settings::computeRotationFuncDesc
bool computeRotationFuncDesc
If true, the rotation function descriptor will be computed, otherwise all its computations will be om...
Definition: ProSHADE_settings.hpp:161
ProSHADE_settings::setSameBoundaries
void setSameBoundaries(bool sameB)
Sets whether same boundaries should be used in the appropriate variable.
Definition: ProSHADE.cpp:619
ProSHADE_settings::maxSymmetryFold
proshade_unsign maxSymmetryFold
The highest symmetry fold to search for.
Definition: ProSHADE_settings.hpp:181
ProSHADE_settings::taylorSeriesCap
proshade_unsign taylorSeriesCap
The max limit on the Taylor series expansion done for the abscissas of the Gauss-Legendre integration...
Definition: ProSHADE_settings.hpp:120
ProSHADE_settings::boundsExtraSpace
proshade_single boundsExtraSpace
The number of extra angstroms to be added to all re-boxing bounds just for safety.
Definition: ProSHADE_settings.hpp:140
ProSHADE_settings::setBicubicInterpolationSearch
void setBicubicInterpolationSearch(bool bicubPeaks)
Sets the bicubic interpolation on peaks.
Definition: ProSHADE.cpp:1148
ProSHADE_settings::setPeakNaiveNoIQR
void setPeakNaiveNoIQR(proshade_double noIQRs)
Sets the number of IQRs from the median for threshold height a peak needs to be considered a peak.
Definition: ProSHADE.cpp:867
ProSHADE_settings::outName
std::string outName
The file name where the output structure(s) should be saved.
Definition: ProSHADE_settings.hpp:155
ProSHADE_settings::setMapInversion
void setMapInversion(bool mInv)
Sets the requested map inversion value in the appropriate variable.
Definition: ProSHADE.cpp:405
ProSHADE_settings::setPeakNeighboursNumber
void setPeakNeighboursNumber(proshade_unsign pkS)
Sets the number of neighbour values that have to be smaller for an index to be considered a peak.
Definition: ProSHADE.cpp:849
ProSHADE_settings::blurFactor
proshade_single blurFactor
This is the amount by which B-factors should be increased to create the blurred map for masking.
Definition: ProSHADE_settings.hpp:129
ProSHADE_settings::maskMap
bool maskMap
Should the map be masked from noise?
Definition: ProSHADE_settings.hpp:131
ProSHADE_settings::requestedSymmetryFold
proshade_unsign requestedSymmetryFold
The fold of the requested symmetry (only applicable to C and D symmetry types).
Definition: ProSHADE_settings.hpp:178
ProSHADE_settings::requestedSymmetryType
std::string requestedSymmetryType
The symmetry type requested by the user. Allowed values are C, D, T, O and I.
Definition: ProSHADE_settings.hpp:177
ProSHADE_settings::setSphereDistances
void setSphereDistances(proshade_single sphDist)
Sets the requested distance between spheres in the appropriate variable.
Definition: ProSHADE.cpp:748
ProSHADE_settings::setVerbosity
void setVerbosity(proshade_signed verbosity)
Sets the requested verbosity in the appropriate variable.
Definition: ProSHADE.cpp:421
ProSHADE_settings::saveMask
bool saveMask
Should the mask be saved?
Definition: ProSHADE_settings.hpp:135
ProSHADE_settings::requestedResolution
proshade_single requestedResolution
The resolution to which the calculations are to be done.
Definition: ProSHADE_settings.hpp:101
ProSHADE_settings::minSymPeak
proshade_double minSymPeak
Minimum average peak for symmetry axis to be considered as "real".
Definition: ProSHADE_settings.hpp:174
ProSHADE_settings::setAxisComparisonThresholdBehaviour
void setAxisComparisonThresholdBehaviour(bool behav)
Sets the automatic symmetry axis tolerance decreasing.
Definition: ProSHADE.cpp:974
ProSHADE_settings::progressiveSphereMapping
bool progressiveSphereMapping
If true, each shell will have its own angular resolution dependent on the actual number of map points...
Definition: ProSHADE_settings.hpp:152
ProSHADE_settings::maxSphereDists
proshade_single maxSphereDists
The distance between spheres in spherical mapping for the largest sphere.
Definition: ProSHADE_settings.hpp:116
ProSHADE_settings::symMissPeakThres
proshade_double symMissPeakThres
Percentage of peaks that could be missing that would warrant starting the missing peaks search proced...
Definition: ProSHADE_settings.hpp:171
ProSHADE_settings::setPDBBFactor
void setPDBBFactor(proshade_double newBF)
Sets the requested B-factor value for PDB files in the appropriate variable.
Definition: ProSHADE.cpp:373
ProSHADE_settings::addStructure
void addStructure(std::string structure)
Adds a structure file name to the appropriate variable.
Definition: ProSHADE.cpp:341
ProSHADE_settings::determineIntegrationOrder
void determineIntegrationOrder(proshade_single maxMapRange)
This function determines the integration order for the between spheres integration.
Definition: ProSHADE.cpp:1275
ProSHADE_settings::setMaskIQR
void setMaskIQR(proshade_single noIQRs)
Sets the requested number of IQRs for masking threshold in the appropriate variable.
Definition: ProSHADE.cpp:454
ProSHADE_settings::setMissingPeakThreshold
void setMissingPeakThreshold(proshade_double mpThres)
Sets the threshold for starting the missing peaks procedure.
Definition: ProSHADE.cpp:938
ProSHADE_settings::maskingThresholdIQRs
proshade_single maskingThresholdIQRs
Number of inter-quartile ranges from the median to be used for thresholding the blurred map for maski...
Definition: ProSHADE_settings.hpp:130
ProSHADE_settings::useCorrelationMasking
bool useCorrelationMasking
Should the blurring masking (false) or the correlation masking (true) be used?
Definition: ProSHADE_settings.hpp:132
ProSHADE_settings::usePhase
bool usePhase
If true, the full data will be used, if false, Patterson maps will be used instead and phased data wi...
Definition: ProSHADE_settings.hpp:113
ProSHADE_settings::setMaskBlurFactor
void setMaskBlurFactor(proshade_single blurFac)
Sets the requested map blurring factor in the appropriate variable.
Definition: ProSHADE.cpp:437
ProSHADE_settings::setPhaseUsage
void setPhaseUsage(bool phaseUsage)
Sets whether the phase information will be used.
Definition: ProSHADE.cpp:885
ProSHADE_settings::maskFileName
std::string maskFileName
The filename to which mask should be saved.
Definition: ProSHADE_settings.hpp:136
ProSHADE_settings::verbose
proshade_signed verbose
Should the software report on the progress, or just be quiet? Value between -1 (nothing) and 4 (loud)
Definition: ProSHADE_settings.hpp:188
ProSHADE_settings::setBandwidth
void setBandwidth(proshade_unsign band)
Sets the requested spherical harmonics bandwidth in the appropriate variable.
Definition: ProSHADE.cpp:732
ProSHADE_settings::setBoundsSpace
void setBoundsSpace(proshade_single boundsExSp)
Sets the requested number of angstroms for extra space in re-boxing in the appropriate variable.
Definition: ProSHADE.cpp:586
ProSHADE_settings::changeMapResolutionTriLinear
bool changeMapResolutionTriLinear
Should maps be re-sampled to obtain the required resolution?
Definition: ProSHADE_settings.hpp:103
ProSHADE_settings::setOutputFilename
void setOutputFilename(std::string oFileName)
Sets the requested output file name in the appropriate variable.
Definition: ProSHADE.cpp:636
ProSHADE_settings::useBiCubicInterpolationOnPeaks
bool useBiCubicInterpolationOnPeaks
This variable switch decides whether best symmetry is detected from peak indices, or whether bicubic ...
Definition: ProSHADE_settings.hpp:180
ProSHADE_internal_spheres::autoDetermineSphereDistances
proshade_single autoDetermineSphereDistances(proshade_single maxMapRange, proshade_single resolution)
This function determines the sphere distances for sphere mapping.
Definition: ProSHADE_spheres.cpp:540
ProSHADE_settings::setAxisComparisonThreshold
void setAxisComparisonThreshold(proshade_double axThres)
Sets the threshold for matching symmetry axes.
Definition: ProSHADE.cpp:955
ProSHADE_settings::forceP1
bool forceP1
Should the P1 spacegroup be forced on the input PDB files?
Definition: ProSHADE_settings.hpp:96
ProSHADE_internal_misc::deepCopyAxisToDblPtrVector
void deepCopyAxisToDblPtrVector(std::vector< proshade_double * > *dblPtrVec, proshade_double *axis)
Does a deep copy of a double array to a vector of double arrays.
Definition: ProSHADE_misc.cpp:287
ProSHADE_settings::task
ProSHADE_Task task
This custom type variable determines which task to perfom (i.e. symmetry detection,...
Definition: ProSHADE_settings.hpp:92
ProSHADE_settings::setMaskFilename
void setMaskFilename(std::string mskFln)
Sets where the mask should be saved.
Definition: ProSHADE.cpp:553
ProSHADE_settings::reBoxMap
bool reBoxMap
This switch decides whether re-boxing is needed.
Definition: ProSHADE_settings.hpp:139
ProSHADE_settings::setSymmetryRotFunPeaks
void setSymmetryRotFunPeaks(bool rotFunPeaks)
Sets the symmetry detection algorithm type.
Definition: ProSHADE.cpp:1134
ProSHADE_settings::normaliseMap
bool normaliseMap
Should the map be normalised to mean 0 sd 1?
Definition: ProSHADE_settings.hpp:123
ProSHADE_settings::addExtraSpace
proshade_single addExtraSpace
If this value is non-zero, this many angstroms of empty space will be added to the internal map.
Definition: ProSHADE_settings.hpp:149
ProSHADE_settings::setRequestedFold
void setRequestedFold(proshade_unsign val)
Sets the user requested symmetry fold.
Definition: ProSHADE.cpp:1060
ProSHADE_settings::determineBandwidthFromAngle
void determineBandwidthFromAngle(proshade_double uncertainty)
This function determines the bandwidth for the spherical harmonics computation from the allowed rotat...
Definition: ProSHADE.cpp:1212
ProSHADE_internal_messages::printWellcomeMessage
void printWellcomeMessage(proshade_signed verbose)
Wellcome message printing.
Definition: ProSHADE_messages.cpp:31
ProSHADE_internal_messages::printTerminateMessage
void printTerminateMessage(proshade_signed verbose)
Final message printing.
Definition: ProSHADE_messages.cpp:49
ProSHADE_settings::setMapReboxing
void setMapReboxing(bool reBx)
Sets whether re-boxing needs to be done in the appropriate variable.
Definition: ProSHADE.cpp:569
ProSHADE_settings::setMaxSymmetryFold
void setMaxSymmetryFold(proshade_unsign maxFold)
Sets the maximum symmetry fold (well, the maximum prime symmetry fold).
Definition: ProSHADE.cpp:1162
ProSHADE_settings::correlationKernel
proshade_single correlationKernel
This value in Angstrom will be used as the kernel for the map-FHM correlation computation.
Definition: ProSHADE_settings.hpp:134
ProSHADE_settings::invertMap
bool invertMap
Should the map be inverted? Only use this if you think you have the wrong hand in your map.
Definition: ProSHADE_settings.hpp:126
ProSHADE_settings::setTaylorSeriesCap
void setTaylorSeriesCap(proshade_unsign tayCap)
Sets the requested Taylor series cap for the Gauss-Legendre integration in the appropriate variable.
Definition: ProSHADE.cpp:781
ProSHADE_settings::boundsSimilarityThreshold
proshade_signed boundsSimilarityThreshold
Number of indices which can be added just to make sure same size in indices is achieved.
Definition: ProSHADE_settings.hpp:141
ProSHADE_settings::setProgressiveSphereMapping
void setProgressiveSphereMapping(bool progSphMap)
Sets the requested sphere mapping value settings approach in the appropriate variable.
Definition: ProSHADE.cpp:716
ProSHADE_settings::changeMapResolution
bool changeMapResolution
Should maps be re-sampled to obtain the required resolution?
Definition: ProSHADE_settings.hpp:102
ProSHADE_settings::setIntegrationOrder
void setIntegrationOrder(proshade_unsign intOrd)
Sets the requested order for the Gauss-Legendre integration in the appropriate variable.
Definition: ProSHADE.cpp:764
ProSHADE_settings::setRotationFunctionComputation
void setRotationFunctionComputation(bool rotfVal)
Sets whether the rotation function distance descriptor should be computed.
Definition: ProSHADE.cpp:832
ProSHADE_settings::rotationUncertainty
proshade_double rotationUncertainty
Alternative to bandwidth - the angle in degrees to which the rotation function accuracy should be com...
Definition: ProSHADE_settings.hpp:110
ProSHADE_settings::setOverlaySaveFile
void setOverlaySaveFile(std::string filename)
Sets the filename to which the overlay structure is to be save into.
Definition: ProSHADE.cpp:1106
ProSHADE_settings::overlayStructureName
std::string overlayStructureName
The filename to which the rotated and translated moving structure is to be saved.
Definition: ProSHADE_settings.hpp:184
ProSHADE_settings::determineBandwidth
void determineBandwidth(proshade_unsign circumference)
This function determines the bandwidth for the spherical harmonics computation.
Definition: ProSHADE.cpp:1180
ProSHADE_settings::moveToCOM
bool moveToCOM
Logical value stating whether the structure should be moved to have its Centre Of Mass (COM) in the m...
Definition: ProSHADE_settings.hpp:146
ProSHADE_settings::peakNeighbours
proshade_unsign peakNeighbours
Number of points in any direction that have to be lower than the considered index in order to conside...
Definition: ProSHADE_settings.hpp:164
ProSHADE_internal_spheres::autoDetermineIntegrationOrder
proshade_unsign autoDetermineIntegrationOrder(proshade_single maxMapRange, proshade_single sphereDist)
This function determines the integration order for the between spheres integration.
Definition: ProSHADE_spheres.cpp:564
ProSHADE_settings::smoothingFactor
proshade_double smoothingFactor
This factor decides how small the group sizes should be - larger factor means more smaller groups.
Definition: ProSHADE_settings.hpp:168
ProSHADE_settings::halfMapKernel
proshade_single halfMapKernel
This value in Angstrom will be used as the kernel for the "fake half-map" computation.
Definition: ProSHADE_settings.hpp:133
ProSHADE_settings::determineSphereDistances
void determineSphereDistances(proshade_single maxMapRange)
This function determines the sphere distances for sphere mapping.
Definition: ProSHADE.cpp:1243
ProSHADE_internal_misc::checkMemoryAllocation
void checkMemoryAllocation(chVar checkVar, std::string fileP, unsigned int lineP, std::string funcP, std::string infoP="This error may occurs when ProSHADE requests memory to be\n : allocated to it and this operation fails. This could\n : happen when not enough memory is available, either due to\n : other processes using a lot of memory, or when the machine\n : does not have sufficient memory available. Re-run to see\n : if this problem persists.")
Checks if memory was allocated properly.
Definition: ProSHADE_misc.hpp:65
ProSHADE_internal_spheres::autoDetermineBandwidth
proshade_unsign autoDetermineBandwidth(proshade_unsign circumference)
This function determines the bandwidth for the spherical harmonics computation.
Definition: ProSHADE_spheres.cpp:518
ProSHADE_settings::setRequestedSymmetry
void setRequestedSymmetry(std::string val)
Sets the user requested symmetry type.
Definition: ProSHADE.cpp:1044
ProSHADE_settings::setExtraSpace
void setExtraSpace(proshade_single exSpace)
Sets the requested map extra space value in the appropriate variable.
Definition: ProSHADE.cpp:700
ProSHADE_settings::forceBounds
proshade_signed * forceBounds
These will be the boundaries to be forced upon the map.
Definition: ProSHADE_settings.hpp:143
ProSHADE_settings::detectedSymmetry
std::vector< proshade_double * > detectedSymmetry
The vector of detected symmetry axes.
Definition: ProSHADE_settings.hpp:192
ProSHADE_settings::recommendedSymmetryFold
proshade_unsign recommendedSymmetryFold
The fold of the recommended symmetry C or D type, 0 otherwise.
Definition: ProSHADE_settings.hpp:176
ProSHADE_settings::inputFiles
std::vector< std::string > inputFiles
This vector contains the filenames of all input structure files.
Definition: ProSHADE_settings.hpp:95
ProSHADE_settings::setMapResolutionChangeTriLinear
void setMapResolutionChangeTriLinear(bool mrChange)
Sets the requested map resolution change decision using tri-linear interpolation in the appropriate v...
Definition: ProSHADE.cpp:668
ProSHADE_settings::setMaskSaving
void setMaskSaving(bool savMsk)
Sets whether the mask should be saved.
Definition: ProSHADE.cpp:537
ProSHADE_settings::firstModelOnly
bool firstModelOnly
Shoud only the first PDB model be used, or should all models be used?
Definition: ProSHADE_settings.hpp:98
ProSHADE_settings::setBoundsThreshold
void setBoundsThreshold(proshade_signed boundsThres)
Sets the threshold for number of indices difference acceptable to make index sizes same in the approp...
Definition: ProSHADE.cpp:602
ProSHADE_settings::axisErrTolerance
proshade_double axisErrTolerance
Allowed error on vector axis in in dot product ( acos ( 1 - axErr ) is the allowed difference in radi...
Definition: ProSHADE_settings.hpp:172
ProSHADE_settings::useSameBounds
bool useSameBounds
Switch to say that the same boundaries as used for the first should be used for all input maps.
Definition: ProSHADE_settings.hpp:142
ProSHADE_settings::setMasking
void setMasking(bool mask)
Sets the requested map masking decision value in the appropriate variable.
Definition: ProSHADE.cpp:470
ProSHADE_settings::pdbBFactorNewVal
proshade_double pdbBFactorNewVal
Change all PDB B-factors to this value (for smooth maps).
Definition: ProSHADE_settings.hpp:106
ProSHADE_settings::usePeakSearchInRotationFunctionSpace
bool usePeakSearchInRotationFunctionSpace
This variable switch decides whether symmetry detection will be done using peak search in rotation fu...
Definition: ProSHADE_settings.hpp:179
ProSHADE_settings::removeWaters
bool removeWaters
Should all waters be removed from input PDB files?
Definition: ProSHADE_settings.hpp:97
ProSHADE_internal_messages::printProgressMessage
void printProgressMessage(proshade_signed verbose, proshade_signed messageLevel, std::string message)
General stdout message printing.
Definition: ProSHADE_messages.cpp:70
ProSHADE_internal_messages::printHelp
void printHelp(void)
This function prints the help screen in the case -h is called, or if command line arguments cannot be...
Definition: ProSHADE_messages.cpp:118
ProSHADE_settings::setResolution
void setResolution(proshade_single resolution)
Sets the requested resolution in the appropriate variable.
Definition: ProSHADE.cpp:357
ProSHADE_internal_misc::addToStringVector
void addToStringVector(std::vector< std::string > *vecToAddTo, std::string elementToAdd)
Adds the element to the vector.
Definition: ProSHADE_misc.cpp:33
ProSHADE_settings::enLevMatrixPowerWeight
proshade_double enLevMatrixPowerWeight
If RRP matrices shell position is to be weighted by putting the position as an exponent,...
Definition: ProSHADE_settings.hpp:159
ProSHADE_settings::setNormalisation
void setNormalisation(bool normalise)
Sets the requested map normalisation value in the appropriate variable.
Definition: ProSHADE.cpp:389
ProSHADE_settings::setEnergyLevelsComputation
void setEnergyLevelsComputation(bool enLevDesc)
Sets whether the energy level distance descriptor should be computed.
Definition: ProSHADE.cpp:798
ProSHADE_settings::computeEnergyLevelsDesc
bool computeEnergyLevelsDesc
If true, the energy levels descriptor will be computed, otherwise all its computations will be omitte...
Definition: ProSHADE_settings.hpp:158