ProSHADE  0.7.6.1 (AUG 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_single xDimAngs, proshade_single yDimAngs, proshade_single 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 setAppliedMaskFilename (std::string mskFln)
 Sets the filename of the mask data that should be applied to the input map. More...
 
void setFourierWeightsFilename (std::string fWgFln)
 Sets the filename of the mask data that should be applied to the input map. 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 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 setFSCThreshold (proshade_double fscThr)
 Sets the minimum FSC threshold for axis to be considered detected. More...
 
void setPeakThreshold (proshade_double peakThr)
 Sets the minimum peak height threshold for axis to be considered possible. More...
 
void setNegativeDensity (bool nDens)
 Sets the internal variable deciding whether input files negative density should be removed. 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?
 
bool removeNegativeDensity
 Should the negative density be removed from input files?
 
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.
 
std::string appliedMaskFileName
 The filename from which mask data will be read from.
 
std::string fourierWeightsFileName
 The filename from which Fourier weights data will be read from.
 
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 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.
 
proshade_double fscThreshold
 The threshold for FSC value under which the axis is considered to be likely noise.
 
proshade_double peakThresholdMin
 The threshold for peak height above which axes are considered possible.
 
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 36 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 37 of file ProSHADE.cpp.

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

◆ 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 160 of file ProSHADE.cpp.

162 {
163  //================================================ Settings regarding the task at hand
164  this->task = taskToPerform;
165 
166  //================================================ Settings regarding input files
167  this->forceP1 = true;
168  this->removeWaters = true;
169  this->firstModelOnly = true;
170  this->removeNegativeDensity = true;
171 
172  //================================================ Settings regarding the resolution of calculations
173  this->requestedResolution = -1.0;
174  this->changeMapResolution = false;
175  this->changeMapResolutionTriLinear = false;
176 
177  //================================================ Settings regarding the PDB B-factor change
178  this->pdbBFactorNewVal = -1.0;
179 
180  //================================================ Settings regarding the bandwidth of calculations
181  this->maxBandwidth = 0;
182  this->rotationUncertainty = 0;
183 
184  //================================================ Settings regarding the phase
185  this->usePhase = true;
186 
187  //================================================ Settings regarding the spheres
188  this->maxSphereDists = 0.0;
189 
190  //================================================ Settings regarding the Gauss-Legendre integration
191  this->integOrder = 0;
192  this->taylorSeriesCap = 10;
193 
194  //================================================ Settings regarding map normalisation
195  this->normaliseMap = false;
196 
197  //================================================ Settings regarding map inversion
198  this->invertMap = false;
199 
200  //================================================ Settings regarding map masking
201  this->blurFactor = 350.0;
202  this->maskingThresholdIQRs = 3.0;
203  this->maskMap = false;
204  this->useCorrelationMasking = false;
205  this->halfMapKernel = 0.0;
206  this->correlationKernel = 0.0;
207  this->saveMask = false;
208  this->maskFileName = "maskFile";
209  this->appliedMaskFileName = "";
210 
211  //================================================ Settings regarding Fourier weights
212  this->fourierWeightsFileName = "";
213 
214  //================================================ Settings regarding re-boxing
215  this->reBoxMap = false;
216  this->boundsExtraSpace = 3.0;
217  this->boundsSimilarityThreshold = 0;
218  this->useSameBounds = false;
219  this->forceBounds = new proshade_signed [6];
220 
221  //================================================ Settings regarding extra cell space
222  this->addExtraSpace = 10.0;
223 
224  //================================================ Settings regarding shell settings
225  this->progressiveSphereMapping = false;
226 
227  //================================================ Settings regarding output file name
228  this->outName = "reBoxed";
229 
230  //================================================ Settings regarding distances computation
231  this->computeEnergyLevelsDesc = true;
232  this->computeTraceSigmaDesc = true;
233  this->computeRotationFuncDesc = true;
234  this->enLevMatrixPowerWeight = 1.0;
235 
236  //================================================ Settings regarding peak searching
237  this->peakNeighbours = 1;
238  this->noIQRsFromMedianNaivePeak = -999.9;
239 
240  //================================================ Settings regarding 1D grouping
241  this->smoothingFactor = 15.0;
242 
243  //================================================ Settings regarding the symmetry detection
244  this->useBiCubicInterpolationOnPeaks = true;
245  this->maxSymmetryFold = 30;
246  this->fscThreshold = 0.65;
247  this->peakThresholdMin = 0.80;
248  this->symMissPeakThres = 0.3;
249  this->axisErrTolerance = 0.01;
250  this->axisErrToleranceDefault = true;
251  this->minSymPeak = 0.3;
252  this->recommendedSymmetryType = "";
253  this->recommendedSymmetryFold = 0;
254  this->requestedSymmetryType = "";
255  this->requestedSymmetryFold = 0;
256  this->detectedSymmetry.clear ( );
257 
258  //================================================ Settings regarding the structure overlay
259  this->overlayStructureName = "movedStructure";
260  this->rotTrsJSONFile = "movedStructureOperations.json";
261 
262  //================================================ Settings regarding verbosity of the program
263  this->verbose = 1;
264 
265  //================================================ Task specific settings
266  switch ( this->task )
267  {
268  case NA:
269  std::cerr << std::endl << "=====================" << std::endl << "!! ProSHADE ERROR !!" << std::endl << "=====================" << std::endl << std::flush;
270  std::cerr << "Error Code : " << "E000014" << std::endl << std::flush;
271  std::cerr << "ProSHADE version : " << PROSHADE_VERSION << std::endl << std::flush;
272  std::cerr << "File : " << "ProSHADE.cpp" << std::endl << std::flush;
273  std::cerr << "Line : " << 97 << std::endl << std::flush;
274  std::cerr << "Function : " << "ProSHADE_settings (Task) constructor" << std::endl << std::flush;
275  std::cerr << "Message : " << "No task has been specified for task specific constructor." << std::endl << std::flush;
276  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;
278  exit ( EXIT_FAILURE );
279 
280  case Symmetry:
281  this->requestedResolution = 6.0;
282  this->pdbBFactorNewVal = 80.0;
283  this->changeMapResolution = true;
284  this->maskMap = false;
285  this->moveToCOM = true;
286  this->reBoxMap = false;
287  break;
288 
289  case Distances:
290  this->changeMapResolution = false;
291  this->maskMap = false;
292  this->moveToCOM = true;
293  this->reBoxMap = false;
294  break;
295 
296  case OverlayMap:
297  this->requestedResolution = 8.0;
298  this->changeMapResolution = true;
299  this->maskMap = false;
300  this->moveToCOM = false;
301  this->reBoxMap = false;
302  break;
303 
304  case MapManip:
305  this->changeMapResolution = false;
306  this->maskMap = true;
307  this->moveToCOM = false;
308  break;
309  }
310 
311  //================================================ Done
312 
313 }

◆ ~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 322 of file ProSHADE.cpp.

324 {
325  //================================================ Release boundaries variable
326  delete[] this->forceBounds;
327 
328  //================================================ Release symmetry axes
329  if ( this->detectedSymmetry.size() > 0 ) { for ( proshade_unsign it = 0; it < static_cast<proshade_unsign> ( this->detectedSymmetry.size() ); it++ ) { if ( this->detectedSymmetry.at(it) != nullptr ) { delete[] this->detectedSymmetry.at(it); } } }
330 
331  //================================================ Done
332 
333 }

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 363 of file ProSHADE.cpp.

365 {
366  //================================================ Use C++ version independent vector processing
367  ProSHADE_internal_misc::addToStringVector ( &( this->inputFiles ), structure );
368 
369  //================================================ Done
370  return ;
371 
372 }

◆ determineAllSHValues()

void ProSHADE_settings::determineAllSHValues ( proshade_unsign  xDim,
proshade_unsign  yDim,
proshade_single  xDimAngs,
proshade_single  yDimAngs,
proshade_single  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]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 1615 of file ProSHADE.cpp.

1616 {
1617  //================================================ Print progress message
1618  ProSHADE_internal_messages::printProgressMessage ( this->verbose, 1, "Preparing spherical harmonics environment." );
1619 
1620  //================================================ Modify dims by resolution
1621  proshade_unsign theoXDim = static_cast< proshade_unsign > ( std::ceil ( xDimAngs / ( this->requestedResolution / 2.0f ) ) );
1622  proshade_unsign theoYDim = static_cast< proshade_unsign > ( std::ceil ( yDimAngs / ( this->requestedResolution / 2.0f ) ) );
1623  proshade_unsign theoZDim = static_cast< proshade_unsign > ( std::ceil ( zDimAngs / ( this->requestedResolution / 2.0f ) ) );
1624 
1625  //================================================ Find maximum circumference
1626  proshade_unsign maxDim = std::max ( theoXDim, std::max ( theoYDim, theoZDim ) );
1627  proshade_unsign minDim = std::min ( theoXDim, std::min ( theoYDim, theoZDim ) );
1628  proshade_unsign midDim = 0;
1629  if ( ( xDim < maxDim ) && ( xDim > minDim ) ) { midDim = theoXDim; }
1630  else if ( ( yDim < maxDim ) && ( yDim > minDim ) ) { midDim = theoYDim; }
1631  else { midDim = theoZDim; }
1632 
1633  proshade_unsign circ = ( maxDim ) + ( midDim );
1634 
1635  //================================================ Bandwidth
1636  if ( this->rotationUncertainty > 0.0 ) { this->determineBandwidthFromAngle ( this->rotationUncertainty ); }
1637  else { this->determineBandwidth ( circ ); }
1638 
1639  //================================================ Find maximum diagonal in Angstroms
1640  proshade_single maxDiag = static_cast< proshade_single > ( std::sqrt ( std::pow ( static_cast<proshade_single> ( maxDim ) * ( this->requestedResolution / 2.0f ), 2.0f ) +
1641  std::pow ( static_cast<proshade_single> ( midDim ) * ( this->requestedResolution / 2.0f ), 2.0f ) ) );
1642 
1643  //================================================ Sphere distances
1644  this->determineSphereDistances ( maxDiag );
1645 
1646  //================================================ Integration order
1647  this->determineIntegrationOrder ( maxDiag );
1648 
1649  //================================================ Report function completion
1650  ProSHADE_internal_messages::printProgressMessage ( this->verbose, 2, "Spherical harmonics environment prepared." );
1651 
1652  //================================================ Done
1653  return ;
1654 
1655 }

◆ 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 1476 of file ProSHADE.cpp.

1477 {
1478  //================================================ Check the current settings value is set to auto
1479  if ( this->maxBandwidth != 0 )
1480  {
1481  std::stringstream hlpSS;
1482  hlpSS << "The bandwidth was determined as: " << this->maxBandwidth;
1484  return ;
1485  }
1486 
1487  //================================================ Determine automatically
1489 
1490  //================================================ Report progress
1491  std::stringstream hlpSS;
1492  hlpSS << "The bandwidth was determined as: " << this->maxBandwidth;
1494 
1495  //================================================ Done
1496  return ;
1497 
1498 }

◆ 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 1508 of file ProSHADE.cpp.

1509 {
1510  //================================================ Determine bandwidth
1511  if ( static_cast<proshade_unsign> ( std::ceil ( ( 360.0 / uncertainty ) / 2 ) ) % 2 == 0 )
1512  {
1513  this->maxBandwidth = static_cast<proshade_unsign> ( std::ceil ( ( 360.0 / uncertainty ) / 2.0 ) );
1514  }
1515  else
1516  {
1517  this->maxBandwidth = static_cast<proshade_unsign> ( std::ceil ( ( 360.0 / uncertainty ) / 2.0 ) ) + 1;
1518  }
1519 
1520  //================================================ Report progress
1521  std::stringstream hlpSS;
1522  hlpSS << "The bandwidth was determined from uncertainty " << uncertainty << " degrees as: " << this->maxBandwidth;
1524 
1525  //================================================ Done
1526  return ;
1527 
1528 }

◆ 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 1571 of file ProSHADE.cpp.

1572 {
1573  //================================================ Check the current settings value is set to auto
1574  if ( this->integOrder != 0 )
1575  {
1576  std::stringstream hlpSS;
1577  hlpSS << "The integration order was determined as " << this->integOrder;
1579  return ;
1580  }
1581 
1582  //================================================ Determine automatically
1584 
1585  //================================================ Report progress
1586  std::stringstream hlpSS;
1587  hlpSS << "The integration order was determined as " << this->integOrder;
1589 
1590  //================================================ Done
1591  return ;
1592 
1593 }

◆ 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 1539 of file ProSHADE.cpp.

1540 {
1541  //================================================ Check the current settings value is set to auto
1542  if ( this->maxSphereDists != 0.0f )
1543  {
1544  std::stringstream hlpSS;
1545  hlpSS << "The sphere distances were determined as " << this->maxSphereDists << " Angstroms.";
1547  return ;
1548  }
1549 
1550  //================================================ Determine automatically
1552 
1553  //================================================ Report progress
1554  std::stringstream hlpSS;
1555  hlpSS << "The sphere distances were determined as " << this->maxSphereDists << " Angstroms.";
1557 
1558  //================================================ Done
1559  return ;
1560 
1561 }

◆ 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 1897 of file ProSHADE.cpp.

1899 {
1900  //================================================ If no command line arguments, print help
1901  if ( argc == 1 ) { ProSHADE_internal_messages::printHelp ( ); }
1902 
1903  //================================================ Long options struct
1904  const struct option_port longopts[] =
1905  {
1906  { "version", no_argument, nullptr, 'v' },
1907  { "help", no_argument, nullptr, 'h' },
1908  { "verbose", required_argument, nullptr, '!' },
1909  { "distances", no_argument, nullptr, 'D' },
1910  { "mapManip", no_argument, nullptr, 'M' },
1911  { "symmetry", no_argument, nullptr, 'S' },
1912  { "overlay", no_argument, nullptr, 'O' },
1913  { "file", required_argument, nullptr, 'f' },
1914  { "forceSpgP1", no_argument, nullptr, 'u' },
1915  { "removeWaters", no_argument, nullptr, 'w' },
1916  { "firstModel", no_argument, nullptr, 'x' },
1917  { "resolution", required_argument, nullptr, 'r' },
1918  { "bandwidth", required_argument, nullptr, 'b' },
1919  { "sphereDists", required_argument, nullptr, 's' },
1920  { "extraSpace", required_argument, nullptr, 'e' },
1921  { "integOrder", required_argument, nullptr, 'i' },
1922  { "taylorCap", required_argument, nullptr, 't' },
1923  { "invertMap", no_argument, nullptr, '@' },
1924  { "normalise", no_argument, nullptr, '#' },
1925  { "mask", no_argument, nullptr, '$' },
1926  { "saveMask", no_argument, nullptr, '%' },
1927  { "maskFile", required_argument, nullptr, '^' },
1928  { "applyMask", required_argument, nullptr, 'G' },
1929  { "maskBlurring", required_argument, nullptr, '&' },
1930  { "maskThreshold", required_argument, nullptr, '*' },
1931  { "mapReboxing", no_argument, nullptr, 'R' },
1932  { "boundsSpace", required_argument, nullptr, '(' },
1933  { "boundsThreshold", required_argument, nullptr, ')' },
1934  { "sameBoundaries", no_argument, nullptr, '-' },
1935  { "reBoxedFilename", required_argument, nullptr, 'g' },
1936  { "pdbTempFact", required_argument, nullptr, 'd' },
1937  { "center", no_argument, nullptr, 'c' },
1938  { "changeMapResol", no_argument, nullptr, 'j' },
1939  { "changeMapTriLin", no_argument, nullptr, 'a' },
1940  { "noPhase", no_argument, nullptr, 'p' },
1941  { "progressive", no_argument, nullptr, 'k' },
1942  { "noEnL", no_argument, nullptr, 'l' },
1943  { "noTrS", no_argument, nullptr, 'm' },
1944  { "noFRF", no_argument, nullptr, 'n' },
1945  { "EnLWeight", required_argument, nullptr, '_' },
1946  { "peakNeigh", required_argument, nullptr, '=' },
1947  { "peakThres", required_argument, nullptr, '+' },
1948  { "missAxThres", required_argument, nullptr, '[' },
1949  { "sameAxComp", required_argument, nullptr, ']' },
1950  { "axisComBeh", no_argument, nullptr, 'q' },
1951  { "bicubSearch", no_argument, nullptr, 'A' },
1952  { "maxSymPrime", required_argument, nullptr, 'B' },
1953  { "minPeakHeight", required_argument, nullptr, 'o' },
1954  { "fscThres", required_argument, nullptr, 'C' },
1955  { "peakMinThres", required_argument, nullptr, 'E' },
1956  { "reqSym", required_argument, nullptr, '{' },
1957  { "overlayFile", required_argument, nullptr, '}' },
1958  { "overlayJSONFile", required_argument, nullptr, 'y' },
1959  { "angUncertain", required_argument, nullptr, ';' },
1960  { "fourierWeights", required_argument, nullptr, 'z' },
1961  { "keepNegDens", no_argument, nullptr, 'F' },
1962  { nullptr, 0, nullptr, 0 }
1963  };
1964 
1965  //================================================ Short options string
1966  const char* const shortopts = "AaB:b:C:cd:DE:e:Ff:G:g:hi:jklmMno:Opqr:Rs:St:uvwxy:z:!:@#$%^:&:*:(:):-_:=:+:[:]:{:}:;:";
1967 
1968  //================================================ Parsing the options
1969  while ( true )
1970  {
1971  //============================================ Read the next option
1972  int opt = getopt_long_port ( argc, argv, shortopts, longopts, nullptr );
1973 
1974  //============================================ Done parsing
1975  if ( opt == -1 )
1976  {
1977  break;
1978  }
1979 
1980  //============================================ For each option, set the internal values appropriately
1981  switch ( opt )
1982  {
1983  //======================================= Print version info
1984  case 'v':
1985  {
1987  exit ( EXIT_SUCCESS );
1988  }
1989 
1990  //======================================= User needs help
1991  case 'h':
1992  {
1994  }
1995 
1996  //======================================= Save the argument as the verbosity value, or if no value given, just set to 3
1997  case '!':
1998  {
1999  this->setVerbosity ( static_cast<proshade_signed> ( atoi ( optarg ) ) );
2000  continue;
2001  }
2002 
2003  //======================================= Set task to distances
2004  case 'D':
2005  {
2006  this->task = Distances;
2007  continue;
2008  }
2009 
2010  //======================================= Set task to map manipulation
2011  case 'M':
2012  {
2013  this->task = MapManip;
2014  continue;
2015  }
2016 
2017  //======================================= Set task to symmetry detection
2018  case 'S':
2019  {
2020  this->task = Symmetry;
2021 
2022  //=================================== Force default unless changed already by the user
2023  const FloatingPoint< proshade_single > lhs1 ( this->requestedResolution ), rhs1 ( -1.0f );
2024  if ( lhs1.AlmostEquals ( rhs1 ) ) { this->requestedResolution = 6.0; }
2025 
2026  const FloatingPoint< proshade_double > lhs2 ( this->pdbBFactorNewVal ), rhs2 ( -1.0 );
2027  if ( lhs2.AlmostEquals ( rhs2 ) ) { this->pdbBFactorNewVal = 80.0; }
2028  this->changeMapResolution = !this->changeMapResolution; // Switch value. This can be over-ridden by the user by using -j
2029  this->moveToCOM = !this->moveToCOM; // Switch value. This can be over-ridden by the user by using -c.
2030 
2031  continue;
2032  }
2033 
2034  //======================================= Set task to map overlay
2035  case 'O':
2036  {
2037  this->task = OverlayMap;
2038  continue;
2039  }
2040 
2041  //======================================= Save the argument as a file to read in
2042  case 'f':
2043  {
2044  this->addStructure ( std::string ( optarg ) );
2045  continue;
2046  }
2047 
2048  //======================================= Force the input PDB files to have P1 spacegroup
2049  case 'u':
2050  {
2051  this->forceP1 = !this->forceP1;
2052  continue;
2053  }
2054 
2055  //======================================= Remove waters from PDB input files?
2056  case 'w':
2057  {
2058  this->removeWaters = !this->removeWaters;
2059  continue;
2060  }
2061 
2062  //======================================= Use all models, or just the first one?
2063  case 'x':
2064  {
2065  this->firstModelOnly = !this->firstModelOnly;
2066  continue;
2067  }
2068 
2069  //======================================= Save the argument as the resolution value
2070  case 'r':
2071  {
2072  this->setResolution ( static_cast<proshade_single> ( atof ( optarg ) ) );
2073  continue;
2074  }
2075 
2076  //======================================= Save the argument as the bandwidth value
2077  case 'b':
2078  {
2079  this->setBandwidth ( static_cast<proshade_unsign> ( atoi ( optarg ) ) );
2080  continue;
2081  }
2082 
2083  //======================================= Save the argument as the extra space value
2084  case 'e':
2085  {
2086  this->setExtraSpace ( static_cast<proshade_single> ( atof ( optarg ) ) );
2087  continue;
2088  }
2089 
2090  //======================================= Save the argument as the intaggration order value
2091  case 'i':
2092  {
2093  this->setIntegrationOrder ( static_cast<proshade_unsign> ( atof ( optarg ) ) );
2094  continue;
2095  }
2096 
2097  //======================================= Save the argument as the sphere distance value
2098  case 's':
2099  {
2100  this->setSphereDistances ( static_cast<proshade_single> ( atof ( optarg ) ) );
2101  continue;
2102  }
2103 
2104  //======================================= Save the argument as the taylor series cap value
2105  case 't':
2106  {
2107  this->setTaylorSeriesCap ( static_cast<proshade_unsign> ( atof ( optarg ) ) );
2108  continue;
2109  }
2110 
2111  //======================================= Set map inversion to true
2112  case '@':
2113  {
2114  this->setMapInversion ( true );
2115  continue;
2116  }
2117 
2118  //======================================= Set map normalisation to true
2119  case '#':
2120  {
2121  this->setNormalisation ( true );
2122  continue;
2123  }
2124 
2125  //======================================= Set map masking to true
2126  case '$':
2127  {
2128  this->setMasking ( true );
2129  continue;
2130  }
2131 
2132  //======================================= Set map masking to true and mask map saving to true as well
2133  case '%':
2134  {
2135  this->setMasking ( true );
2136  this->setMaskSaving ( true );
2137  continue;
2138  }
2139 
2140  //======================================= Save the argument as the mask filename value
2141  case '^':
2142  {
2143  this->setMaskFilename ( static_cast<std::string> ( optarg ) );
2144  continue;
2145  }
2146 
2147  //======================================= Save the argument as the mask filename value
2148  case 'G':
2149  {
2150  this->setAppliedMaskFilename ( static_cast<std::string> ( optarg ) );
2151  continue;
2152  }
2153 
2154  //======================================= Save the argument as the Fourier weights filename value
2155  case 'z':
2156  {
2157  this->setFourierWeightsFilename ( static_cast<std::string> ( optarg ) );
2158  continue;
2159  }
2160 
2161  //======================================= Save the argument as the mask blurring factor value
2162  case '&':
2163  {
2164  this->setMaskBlurFactor ( static_cast<proshade_single> ( atof ( optarg ) ) );
2165  continue;
2166  }
2167 
2168  //======================================= Save the argument as the mask threshold (IQR) value
2169  case '*':
2170  {
2171  this->setMaskIQR ( static_cast<proshade_single> ( atof ( optarg ) ) );
2172  continue;
2173  }
2174 
2175  //======================================= Set map reboxing to true
2176  case 'R':
2177  {
2178  this->setMasking ( true );
2179  this->setMapReboxing ( true );
2180  continue;
2181  }
2182 
2183  //======================================= Save the argument as the bounds extra space value
2184  case '(':
2185  {
2186  this->setBoundsSpace ( static_cast<proshade_single> ( atof ( optarg ) ) );
2187  continue;
2188  }
2189 
2190  //======================================= Save the argument as the bounds similarity threshold value
2191  case ')':
2192  {
2193  this->setBoundsThreshold ( static_cast<proshade_signed> ( atoi ( optarg ) ) );
2194  continue;
2195  }
2196 
2197  //======================================= Set same boundaries to true
2198  case '-':
2199  {
2200  this->setSameBoundaries ( true );
2201  continue;
2202  }
2203 
2204  //======================================= Save the argument as the re-boxed structure filename value
2205  case 'g':
2206  {
2207  this->setOutputFilename ( static_cast<std::string> ( optarg ) );
2208  continue;
2209  }
2210 
2211  //======================================= Save the argument as the PDB B-factor new constant value
2212  case 'd':
2213  {
2214  this->setPDBBFactor ( static_cast<proshade_double> ( atof ( optarg ) ) );
2215  continue;
2216  }
2217 
2218  //======================================= Set map centering to true
2219  case 'c':
2220  {
2221  this->moveToCOM = !this->moveToCOM;
2222  continue;
2223  }
2224 
2225  //======================================= Set map resolution change using Fourier transforms to true
2226  case 'j':
2227  {
2229  continue;
2230  }
2231 
2232  //======================================= Set map resolution change using real-space tri-linear interpolation to true
2233  case 'a':
2234  {
2235  this->setMapResolutionChangeTriLinear ( true );
2236  continue;
2237  }
2238 
2239  //======================================= Set map phase removal to true
2240  case 'p':
2241  {
2242  this->setPhaseUsage ( false );
2243  continue;
2244  }
2245 
2246  //======================================= Set progressive shell mapping to true
2247  case 'k':
2248  {
2249  this->setProgressiveSphereMapping ( true );
2250  continue;
2251  }
2252 
2253  //======================================= Set energy level descriptor computation to false
2254  case 'l':
2255  {
2256  this->setEnergyLevelsComputation ( false );
2257  continue;
2258  }
2259 
2260  //======================================= Set trace sigma descriptor computation to false
2261  case 'm':
2262  {
2263  this->setTraceSigmaComputation ( false );
2264  continue;
2265  }
2266 
2267  //======================================= Set full rotation function descriptor computation to false
2268  case 'n':
2269  {
2270  this->setRotationFunctionComputation ( false );
2271  continue;
2272  }
2273 
2274  //======================================= Save the argument as the energy levels descriptor weight value
2275  case '_':
2276  {
2277  this->setEnLevShellWeight ( static_cast<proshade_double> ( atof ( optarg ) ) );
2278  continue;
2279  }
2280 
2281  //======================================= Save the argument as the peak neighbours minimum value
2282  case '=':
2283  {
2284  this->setPeakNeighboursNumber ( static_cast<proshade_unsign> ( atoi ( optarg ) ) );
2285  continue;
2286  }
2287 
2288  //======================================= Save the argument as the peak IQR from median naive small peaks removal value
2289  case '+':
2290  {
2291  this->setPeakNaiveNoIQR ( static_cast<proshade_double> ( atof ( optarg ) ) );
2292  continue;
2293  }
2294 
2295  //======================================= Save the argument as the missing axis threshold value
2296  case '[':
2297  {
2298  this->setMissingPeakThreshold ( static_cast<proshade_double> ( atof ( optarg ) ) );
2299  continue;
2300  }
2301 
2302  //======================================= Save the argument as the missing axis threshold value
2303  case ']':
2304  {
2305  setAxisComparisonThreshold ( static_cast<proshade_double> ( atof ( optarg ) ) );
2306  continue;
2307  }
2308 
2309  //======================================= Save the argument as the missing axis threshold value
2310  case 'q':
2311  {
2312  this->setAxisComparisonThresholdBehaviour ( !this->axisErrToleranceDefault );
2313  continue;
2314  }
2315 
2316  //======================================= Save the argument as the bicubic interpolation search requirement value
2317  case 'A':
2318  {
2320  continue;
2321  }
2322 
2323  //======================================= Save the argument as the bicubic interpolation search requirement value
2324  case 'B':
2325  {
2326  setMaxSymmetryFold ( static_cast<proshade_unsign> ( atoi ( optarg ) ) );
2327  continue;
2328  }
2329 
2330  //======================================= Minimum peak height for axis
2331  case 'o':
2332  {
2333  this->minSymPeak = static_cast<proshade_double> ( atof ( optarg ) );
2334  continue;
2335  }
2336 
2337  //======================================= Minimum FSC value for axis to be detected
2338  case 'C':
2339  {
2340  this->setFSCThreshold ( static_cast<proshade_double> ( atof ( optarg ) ) );
2341  continue;
2342  }
2343 
2344  //======================================= Minimum peak height value for axis to be considered possible
2345  case 'E':
2346  {
2347  this->setPeakThreshold ( static_cast<proshade_double> ( atof ( optarg ) ) );
2348  continue;
2349  }
2350 
2351  //======================================= Save the argument as the requested symmetry and potentially fold value
2352  case '{':
2353  {
2354  std::string input = static_cast<std::string> ( optarg );
2355 
2356  if ( input.at(0) == 'C' )
2357  {
2358  this->setRequestedSymmetry ( "C" );
2359 
2360  std::string numHlp ( input.begin()+1, input.end() );
2361  if ( numHlp.length() > 0 ) { this->setRequestedFold ( static_cast< proshade_unsign > ( atoi ( numHlp.c_str() ) ) ); }
2362  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 ); }
2363  }
2364  else
2365  {
2366  if ( input.at(0) == 'D' )
2367  {
2368  this->setRequestedSymmetry ( "D" );
2369 
2370  std::string numHlp ( input.begin()+1, input.end() );
2371  if ( numHlp.length() > 0 ) { this->setRequestedFold ( static_cast< proshade_unsign > ( atoi ( numHlp.c_str() ) ) ); }
2372  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 ); }
2373  }
2374  else
2375  {
2376  if ( input.at(0) == 'T' )
2377  {
2378  this->setRequestedSymmetry ( "T" );
2379  }
2380  else
2381  {
2382  if ( input.at(0) == 'O' )
2383  {
2384  this->setRequestedSymmetry ( "O" );
2385  }
2386  else
2387  {
2388  if ( input.at(0) == 'I' )
2389  {
2390  this->setRequestedSymmetry ( "I" );
2391  }
2392  else
2393  {
2394  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 );
2395  }
2396  }
2397  }
2398  }
2399  }
2400 
2401  continue;
2402  }
2403 
2404  //======================================= Save the argument as filename to save the overlay moved structure to value
2405  case '}':
2406  {
2407  this->setOverlaySaveFile ( static_cast<std::string> ( optarg ) );
2408  continue;
2409  }
2410 
2411  //======================================= Save the argument as filename to save the overlay operations to value
2412  case 'y':
2413  {
2414  this->setOverlayJsonFile ( static_cast<std::string> ( optarg ) );
2415  continue;
2416  }
2417 
2418  //======================================= Save the argument as angular uncertainty for bandwidth determination
2419  case ';':
2420  {
2421  this->rotationUncertainty = static_cast<proshade_double> ( atof ( optarg ) );
2422  continue;
2423  }
2424 
2425  //======================================= Should the negative density from input files be removed?
2426  case 'F':
2427  {
2428  this->setNegativeDensity ( false );
2429  continue;
2430  }
2431 
2432  //======================================= Unknown option
2433  case '?':
2434  {
2435  //=================================== Save the argument as angular uncertainty for bandwidth determination
2436  if ( optopt )
2437  {
2438  std::cout << "!!! ProSHADE ERROR !!! Unrecognised short option -" << static_cast<char> ( optopt ) << " . Please use -h for help on the command line options." << std::endl;
2439  }
2440  else
2441  {
2442  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;
2443  }
2444 
2445  //=================================== This case is handled by getopt_long, nothing more needed.
2446  exit ( EXIT_SUCCESS );
2447  }
2448 
2449  //======================================= Fallback option
2450  default:
2451  {
2453  }
2454  }
2455  }
2456 
2457  //================================================ Done
2458  return ;
2459 
2460 }

◆ 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 2469 of file ProSHADE.cpp.

2471 {
2472  //================================================ Print the currest values in the settings object
2473  //== Settings regarding the task at hand
2474  std::stringstream strstr;
2475  strstr.str(std::string());
2476  if ( this->task == NA ) { strstr << "NA"; }
2477  if ( this->task == Distances ) { strstr << "DISTANCES COMPUTATION"; }
2478  if ( this->task == MapManip ) { strstr << "MAP MANIPULATION"; }
2479  if ( this->task == Symmetry ) { strstr << "SYMMETRY DETECTION"; }
2480  if ( this->task == OverlayMap ) { strstr << "MAP OVERLAY"; }
2481  printf ( "Task to perform : %37s\n", strstr.str().c_str() );
2482 
2483  //== Settings regarding the input files
2484  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( this->inputFiles.size() ); iter++ )
2485  {
2486  strstr.str(std::string());
2487  strstr << this->inputFiles.at(iter);
2488  printf ( "File(s) to process : %37s\n", strstr.str().c_str() );
2489  }
2490 
2491  strstr.str(std::string());
2492  if ( this->forceP1 ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2493  printf ( "Force P1 spacegroup : %37s\n", strstr.str().c_str() );
2494 
2495  strstr.str(std::string());
2496  if ( this->removeWaters ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2497  printf ( "Waters removed : %37s\n", strstr.str().c_str() );
2498 
2499  strstr.str(std::string());
2500  if ( this->firstModelOnly ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2501  printf ( "Only 1st model : %37s\n", strstr.str().c_str() );
2502 
2503  strstr.str(std::string());
2504  if ( this->removeNegativeDensity ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2505  printf ( "Remove neg. dens. : %37s\n", strstr.str().c_str() );
2506 
2507  //== Settings regarding the resolution of calculations
2508  strstr.str(std::string());
2509  strstr << this->requestedResolution;
2510  printf ( "Resolution (comp) : %37s\n", strstr.str().c_str() );
2511 
2512  strstr.str(std::string());
2513  if ( this->changeMapResolution ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2514  printf ( "Change map resol : %37s\n", strstr.str().c_str() );
2515 
2516  strstr.str(std::string());
2517  if ( this->changeMapResolutionTriLinear ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2518  printf ( "Change map tri-lin : %37s\n", strstr.str().c_str() );
2519 
2520  //== Settings regarding the PDB B-factor change
2521  strstr.str(std::string());
2522  strstr << this->pdbBFactorNewVal;
2523  printf ( "PDB B-factor const : %37s\n", strstr.str().c_str() );
2524 
2525  //== Settings regarding the bandwidth of calculations
2526  strstr.str(std::string());
2527  strstr << this->maxBandwidth;
2528  printf ( "Bandwidth : %37s\n", strstr.str().c_str() );
2529 
2530  strstr.str(std::string());
2531  strstr << this->rotationUncertainty;
2532  printf ( "Rotation doubt : %37s\n", strstr.str().c_str() );
2533 
2534  //== Settings regarding the phase
2535  strstr.str(std::string());
2536  if ( this->usePhase ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2537  printf ( "Use phase info : %37s\n", strstr.str().c_str() );
2538 
2539  //== Settings regarding the spheres
2540  strstr.str(std::string());
2541  strstr << this->maxSphereDists;
2542  printf ( "Sphere distances : %37s\n", strstr.str().c_str() );
2543 
2544  //== Settings regarding the Gauss-Legendre integration
2545  strstr.str(std::string());
2546  strstr << this->integOrder;
2547  printf ( "Integration order : %37s\n", strstr.str().c_str() );
2548 
2549  strstr.str(std::string());
2550  strstr << this->taylorSeriesCap;
2551  printf ( "Taylor series cap : %37s\n", strstr.str().c_str() );
2552 
2553  //== Settings regarding map normalisation
2554  strstr.str(std::string());
2555  if ( this->normaliseMap ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2556  printf ( "Map normalisation : %37s\n", strstr.str().c_str() );
2557 
2558  //== Settings regarding map inversion
2559  strstr.str(std::string());
2560  if ( this->invertMap ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2561  printf ( "Map inversion : %37s\n", strstr.str().c_str() );
2562 
2563  //== Settings regarding map masking
2564  strstr.str(std::string());
2565  strstr << this->blurFactor;
2566  printf ( "Map blurring : %37s\n", strstr.str().c_str() );
2567 
2568  strstr.str(std::string());
2569  strstr << this->maskingThresholdIQRs;
2570  printf ( "Masking threshold : %37s\n", strstr.str().c_str() );
2571 
2572  strstr.str(std::string());
2573  if ( this->maskMap ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2574  printf ( "Map masking : %37s\n", strstr.str().c_str() );
2575 
2576  strstr.str(std::string());
2577  if ( this->useCorrelationMasking ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2578  printf ( "Correlation mask : %37s\n", strstr.str().c_str() );
2579 
2580  strstr.str(std::string());
2581  strstr << this->halfMapKernel;
2582  printf ( "Half-map kernel : %37s\n", strstr.str().c_str() );
2583 
2584  strstr.str(std::string());
2585  strstr << this->correlationKernel;
2586  printf ( "Correlation kernel : %37s\n", strstr.str().c_str() );
2587 
2588  strstr.str(std::string());
2589  if ( this->saveMask ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2590  printf ( "Saving mask : %37s\n", strstr.str().c_str() );
2591 
2592  strstr.str(std::string());
2593  strstr << this->maskFileName;
2594  printf ( "Map mask filename : %37s\n", strstr.str().c_str() );
2595 
2596  //== Settings regarding re-boxing
2597  strstr.str(std::string());
2598  if ( this->reBoxMap ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2599  printf ( "Map re-boxing : %37s\n", strstr.str().c_str() );
2600 
2601  strstr.str(std::string());
2602  strstr << this->boundsExtraSpace;
2603  printf ( "Bounds extra space : %37s\n", strstr.str().c_str() );
2604 
2605  strstr.str(std::string());
2606  strstr << this->boundsSimilarityThreshold;
2607  printf ( "Bounds similarity : %37s\n", strstr.str().c_str() );
2608 
2609  strstr.str(std::string());
2610  if ( this->useSameBounds ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2611  printf ( "Same boundaries : %37s\n", strstr.str().c_str() );
2612 
2613  strstr.str(std::string());
2614  if ( this->forceBounds != nullptr )
2615  {
2616  strstr << this->forceBounds[0] << " - " << this->forceBounds[1] << " | " << this->forceBounds[2] << " - " << this->forceBounds[3] << " | " << this->forceBounds[4] << " - " << this->forceBounds[5];
2617  printf ( "Bounds similarity : %37s\n", strstr.str().c_str() );
2618  }
2619  else
2620  {
2621  strstr << "Not allocated";
2622  printf ( "Bounds similarity : %37s\n", strstr.str().c_str() );
2623  }
2624 
2625  //== Settings regarding COM
2626  strstr.str(std::string());
2627  if ( this->moveToCOM ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2628  printf ( "Map COM centering : %37s\n", strstr.str().c_str() );
2629 
2630  //== Settings regarding extra cell space
2631  strstr.str(std::string());
2632  strstr << this->addExtraSpace;
2633  printf ( "Extra space : %37s\n", strstr.str().c_str() );
2634 
2635  //== Settings regarding shell settings
2636  strstr.str(std::string());
2637  if ( this->progressiveSphereMapping ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2638  printf ( "Progressive spheres : %37s\n", strstr.str().c_str() );
2639 
2640  //== Settings regarding output file name
2641  strstr.str(std::string());
2642  strstr << this->outName;
2643  printf ( "Re-boxed filename : %37s\n", strstr.str().c_str() );
2644 
2645  //== Settings regarding distances computation
2646  strstr.str(std::string());
2647  if ( this->computeEnergyLevelsDesc ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2648  printf ( "Energy lvl desc : %37s\n", strstr.str().c_str() );
2649 
2650  strstr.str(std::string());
2651  strstr << this->enLevMatrixPowerWeight;
2652  printf ( "Energy lvl weight : %37s\n", strstr.str().c_str() );
2653 
2654  strstr.str(std::string());
2655  if ( this->computeTraceSigmaDesc ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2656  printf ( "Tr sigma desc : %37s\n", strstr.str().c_str() );
2657 
2658  strstr.str(std::string());
2659  if ( this->computeRotationFuncDesc ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2660  printf ( "Full RF desc : %37s\n", strstr.str().c_str() );
2661 
2662  //== Settings regarding peak searching
2663  strstr.str(std::string());
2664  strstr << this->peakNeighbours;
2665  printf ( "Neightbours to peak : %37s\n", strstr.str().c_str() );
2666 
2667  strstr.str(std::string());
2668  strstr << this->noIQRsFromMedianNaivePeak;
2669  printf ( "Peak IQR threshold : %37s\n", strstr.str().c_str() );
2670 
2671  //== Settings regarding 1D grouping
2672  strstr.str(std::string());
2673  strstr << this->smoothingFactor;
2674  printf ( "Smoothing factor : %37s\n", strstr.str().c_str() );
2675 
2676  //== Settings regarding the symmetry detection
2677  strstr.str(std::string());
2678  strstr << this->symMissPeakThres;
2679  printf ( "Missing ax. thres : %37s\n", strstr.str().c_str() );
2680 
2681  strstr.str(std::string());
2682  strstr << this->axisErrTolerance;
2683  printf ( "Same ax. threshold : %37s\n", strstr.str().c_str() );
2684 
2685  strstr.str(std::string());
2686  if ( this->axisErrToleranceDefault ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2687  printf ( "Same ax. thre. decr.: %37s\n", strstr.str().c_str() );
2688 
2689  strstr.str(std::string());
2690  strstr << this->minSymPeak;
2691  printf ( "Min. sym. peak size : %37s\n", strstr.str().c_str() );
2692 
2693  strstr.str(std::string());
2694  strstr << this->recommendedSymmetryType << "-" << this->recommendedSymmetryFold;
2695  printf ( "Recommended symm. : %37s\n", strstr.str().c_str() );
2696 
2697  strstr.str(std::string());
2698  strstr << this->requestedSymmetryType << "-" << this->requestedSymmetryFold;
2699  printf ( "Requested symm. : %37s\n", strstr.str().c_str() );
2700 
2701  strstr.str(std::string());
2702  if ( this->useBiCubicInterpolationOnPeaks ) { strstr << "TRUE"; } else { strstr << "FALSE"; }
2703  printf ( "Use bicubic interp. : %37s\n", strstr.str().c_str() );
2704 
2705  strstr.str(std::string());
2706  strstr << this->maxSymmetryFold;
2707  printf ( "Max symmetry fold : %37s\n", strstr.str().c_str() );
2708 
2709  strstr.str(std::string());
2710  strstr << this->fscThreshold;
2711  printf ( "FSC Threshold : %37s\n", strstr.str().c_str() );
2712 
2713  strstr.str(std::string());
2714  strstr << this->peakThresholdMin;
2715  printf ( "Peak Threshold : %37s\n", strstr.str().c_str() );
2716 
2717  //== Settings regarding the structure overlay
2718  strstr.str(std::string());
2719  strstr << this->overlayStructureName;
2720  printf ( "Overlay file : %37s\n", strstr.str().c_str() );
2721 
2722  strstr.str(std::string());
2723  strstr << this->rotTrsJSONFile;
2724  printf ( "JSON overlay file : %37s\n", strstr.str().c_str() );
2725 
2726  //== Settings regarding verbosity of the program
2727  strstr.str(std::string());
2728  strstr << this->verbose;
2729  printf ( "Verbosity : %37s\n", strstr.str().c_str() );
2730 
2731  //================================================ Done
2732  return ;
2733 
2734 }

◆ setAppliedMaskFilename()

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

Sets the filename of the mask data that should be applied to the input map.

This function sets the the filename from which mask should be read from.

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

Definition at line 647 of file ProSHADE.cpp.

649 {
650  //================================================ Set the value
651  this->appliedMaskFileName = mskFln;
652 
653  //================================================ Done
654  return ;
655 
656 }

◆ 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 1165 of file ProSHADE.cpp.

1167 {
1168  //================================================ Set the value
1169  this->axisErrTolerance = axThres;
1170 
1171  //================================================ Done
1172  return ;
1173 
1174 }

◆ 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 1188 of file ProSHADE.cpp.

1190 {
1191  //================================================ Set the value
1192  this->axisErrToleranceDefault = behav;
1193 
1194  //================================================ Done
1195  return ;
1196 
1197 }

◆ 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 870 of file ProSHADE.cpp.

872 {
873  //================================================ Set the value
874  this->maxBandwidth = band;
875 
876  //================================================ Done
877  return ;
878 
879 }

◆ 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 1385 of file ProSHADE.cpp.

1387 {
1388  //================================================ Set the value
1389  this->useBiCubicInterpolationOnPeaks = bicubPeaks;
1390 
1391  //================================================ Done
1392  return ;
1393 
1394 }

◆ 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 708 of file ProSHADE.cpp.

710 {
711  //================================================ Set the value
712  this->boundsExtraSpace = boundsExSp;
713 
714  //================================================ Done
715  return ;
716 
717 }

◆ 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 728 of file ProSHADE.cpp.

730 {
731  //================================================ Set the value
732  this->boundsSimilarityThreshold = boundsThres;
733 
734  //================================================ Done
735  return ;
736 
737 }

◆ 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 545 of file ProSHADE.cpp.

547 {
548  //================================================ Set the value
549  this->useCorrelationMasking = corMask;
550 
551  //================================================ Done
552  return ;
553 
554 }

◆ 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 1315 of file ProSHADE.cpp.

1317 {
1318  //================================================ Allocate memory
1319  proshade_double* hlpAxis = new proshade_double [7];
1320  ProSHADE_internal_misc::checkMemoryAllocation ( hlpAxis, __FILE__, __LINE__, __func__ );
1321 
1322  //================================================ Copy (deep) data
1323  hlpAxis[0] = sym[0];
1324  hlpAxis[1] = sym[1];
1325  hlpAxis[2] = sym[2];
1326  hlpAxis[3] = sym[3];
1327  hlpAxis[4] = sym[4];
1328  hlpAxis[5] = sym[5];
1329  hlpAxis[6] = sym[6];
1330 
1331  //================================================ Save
1333 
1334  //================================================ Release memory
1335  delete[] hlpAxis;
1336 
1337  //================================================ Done
1338  return ;
1339 
1340 }

◆ 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 972 of file ProSHADE.cpp.

974 {
975  //======================================== Set the value
976  this->computeEnergyLevelsDesc = enLevDesc;
977 
978  //======================================== Done
979  return ;
980 
981 }

◆ 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 1101 of file ProSHADE.cpp.

1103 {
1104  //================================================ Set the value
1105  this->enLevMatrixPowerWeight = mPower;
1106 
1107  //================================================ Done
1108  return ;
1109 
1110 }

◆ 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 850 of file ProSHADE.cpp.

852 {
853  //================================================ Set the value
854  this->addExtraSpace = exSpace;
855 
856  //================================================ Done
857  return ;
858 
859 }

◆ setFourierWeightsFilename()

void ProSHADE_settings::setFourierWeightsFilename ( std::string  fWgFln)

Sets the filename of the mask data that should be applied to the input map.

This function sets the the filename from which mask should be read from.

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

Definition at line 667 of file ProSHADE.cpp.

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

◆ setFSCThreshold()

void ProSHADE_settings::setFSCThreshold ( proshade_double  fscThr)

Sets the minimum FSC threshold for axis to be considered detected.

Parameters
[in]fscThrThe minimum axis FSC threshold for the axis to be considered detected.

Definition at line 1421 of file ProSHADE.cpp.

1423 {
1424  //================================================ Set the value
1425  this->fscThreshold = fscThr;
1426 
1427  //================================================ Done
1428  return ;
1429 
1430 }

◆ 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 1123 of file ProSHADE.cpp.

1125 {
1126  //================================================ Set the value
1127  this->smoothingFactor = smFact;
1128 
1129  //================================================ Done
1130  return ;
1131 
1132 }

◆ 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 910 of file ProSHADE.cpp.

912 {
913  //================================================ Set the value
914  this->integOrder = intOrd;
915 
916  //================================================ Done
917  return ;
918 
919 }

◆ 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 830 of file ProSHADE.cpp.

832 {
833  //================================================ Set the value
834  this->moveToCOM = com;
835 
836  //================================================ Done
837  return ;
838 
839 }

◆ 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 443 of file ProSHADE.cpp.

445 {
446  //================================================ Set the value
447  this->invertMap = mInv;
448 
449  //================================================ Done
450  return ;
451 
452 }

◆ 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 687 of file ProSHADE.cpp.

689 {
690  //================================================ Set the value
691  this->reBoxMap = reBx;
692 
693  //================================================ Done
694  return ;
695 
696 }

◆ 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 790 of file ProSHADE.cpp.

792 {
793  //================================================ Set the value
794  this->changeMapResolution = mrChange;
795 
796  //================================================ Done
797  return ;
798 
799 }

◆ 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 810 of file ProSHADE.cpp.

812 {
813  //================================================ Set the value
814  this->changeMapResolutionTriLinear = mrChange;
815 
816  //================================================ Done
817  return ;
818 
819 }

◆ 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 483 of file ProSHADE.cpp.

485 {
486  //================================================ Set the value
487  this->blurFactor = blurFac;
488 
489  //================================================ Done
490  return ;
491 
492 }

◆ 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 627 of file ProSHADE.cpp.

629 {
630  //================================================ Set the value
631  this->maskFileName = mskFln;
632 
633  //================================================ Done
634  return ;
635 
636 }

◆ 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 524 of file ProSHADE.cpp.

526 {
527  //================================================ Set the value
528  this->maskMap = mask;
529 
530  //================================================ Done
531  return ;
532 
533 }

◆ 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 504 of file ProSHADE.cpp.

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

◆ 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 607 of file ProSHADE.cpp.

609 {
610  //================================================ Set the value
611  this->saveMask = savMsk;
612 
613  //================================================ Done
614  return ;
615 
616 }

◆ 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 1403 of file ProSHADE.cpp.

1405 {
1406  //================================================ Set the value
1407  this->maxSymmetryFold = maxFold;
1408 
1409  //================================================ Done
1410  return ;
1411 
1412 }

◆ 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 587 of file ProSHADE.cpp.

589 {
590  //================================================ Set the value
591  this->correlationKernel = minMS;
592 
593  //================================================ Done
594  return ;
595 
596 }

◆ 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 1209 of file ProSHADE.cpp.

1211 {
1212  //================================================ Set the value
1213  this->minSymPeak = minSP;
1214 
1215  //================================================ Done
1216  return ;
1217 
1218 }

◆ 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 1144 of file ProSHADE.cpp.

1146 {
1147  //================================================ Set the value
1148  this->symMissPeakThres = mpThres;
1149 
1150  //================================================ Done
1151  return ;
1152 
1153 }

◆ setNegativeDensity()

void ProSHADE_settings::setNegativeDensity ( bool  nDens)

Sets the internal variable deciding whether input files negative density should be removed.

Parameters
[in]nDensShould the negative density be removed from input files?

Definition at line 1457 of file ProSHADE.cpp.

1459 {
1460  //================================================ Set the value
1461  this->removeNegativeDensity = nDens;
1462 
1463  //================================================ Done
1464  return ;
1465 
1466 }

◆ 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 423 of file ProSHADE.cpp.

425 {
426  //================================================ Set the value
427  this->normaliseMap = normalise;
428 
429  //================================================ Done
430  return ;
431 
432 }

◆ 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 770 of file ProSHADE.cpp.

772 {
773  //================================================ Set the value
774  this->outName = oFileName;
775 
776  //================================================ Done
777  return ;
778 
779 }

◆ 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 1367 of file ProSHADE.cpp.

1369 {
1370  //================================================ Set the value
1371  this->rotTrsJSONFile = filename;
1372 
1373  //================================================ Done
1374  return ;
1375 
1376 }

◆ 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 1349 of file ProSHADE.cpp.

1351 {
1352  //================================================ Set the value
1353  this->overlayStructureName = filename;
1354 
1355  //================================================ Done
1356  return ;
1357 
1358 }

◆ 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 403 of file ProSHADE.cpp.

405 {
406  //================================================ Set the value
407  this->pdbBFactorNewVal = newBF;
408 
409  //================================================ Done
410  return ;
411 
412 }

◆ 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 1057 of file ProSHADE.cpp.

1059 {
1060  //================================================ Set the value
1061  this->noIQRsFromMedianNaivePeak = noIQRs;
1062 
1063  //================================================ Done
1064  return ;
1065 
1066 }

◆ 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 1035 of file ProSHADE.cpp.

1037 {
1038  //================================================ Set the value
1039  this->peakNeighbours = pkS;
1040 
1041  //================================================ Done
1042  return ;
1043 
1044 }

◆ setPeakThreshold()

void ProSHADE_settings::setPeakThreshold ( proshade_double  peakThr)

Sets the minimum peak height threshold for axis to be considered possible.

Parameters
[in]fscThrThe minimum axis peak height threshold for the axis to be considered possible.

Definition at line 1439 of file ProSHADE.cpp.

1441 {
1442  //================================================ Set the value
1443  this->peakThresholdMin = peakThr;
1444 
1445  //================================================ Done
1446  return ;
1447 
1448 }

◆ 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 1079 of file ProSHADE.cpp.

1081 {
1082  //================================================ Set the value
1083  this->usePhase = phaseUsage;
1084 
1085  //================================================ Done
1086  return ;
1087 
1088 }

◆ 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 951 of file ProSHADE.cpp.

953 {
954  //================================================ Set the value
955  this->progressiveSphereMapping = progSphMap;
956 
957  //================================================ Done
958  return ;
959 
960 }

◆ 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 1254 of file ProSHADE.cpp.

1256 {
1257  //================================================ Set the value
1258  this->recommendedSymmetryFold = val;
1259 
1260  //================================================ Done
1261  return ;
1262 
1263 }

◆ 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 1231 of file ProSHADE.cpp.

1233 {
1234  //================================================ Set the value
1235  this->recommendedSymmetryType = val;
1236 
1237  //================================================ Done
1238  return ;
1239 
1240 }

◆ 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 1294 of file ProSHADE.cpp.

1296 {
1297  //================================================ Set the value
1298  this->requestedSymmetryFold = val;
1299 
1300  //================================================ Done
1301  return ;
1302 
1303 }

◆ 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 1274 of file ProSHADE.cpp.

1276 {
1277  //================================================ Set the value
1278  this->requestedSymmetryType = val;
1279 
1280  //================================================ Done
1281  return ;
1282 
1283 }

◆ 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 383 of file ProSHADE.cpp.

385 {
386  //================================================ Set the value
387  this->requestedResolution = resolution;
388 
389  //================================================ Done
390  return ;
391 
392 }

◆ 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 1014 of file ProSHADE.cpp.

1016 {
1017  //================================================ Set the value
1018  this->computeRotationFuncDesc = rotfVal;
1019 
1020  //================================================ Done
1021  return ;
1022 
1023 }

◆ 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 749 of file ProSHADE.cpp.

751 {
752  //================================================ Set the value
753  this->useSameBounds = sameB;
754 
755  //================================================ Done
756  return ;
757 
758 }

◆ 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 890 of file ProSHADE.cpp.

892 {
893  //================================================ Set the value
894  this->maxSphereDists = sphDist;
895 
896  //================================================ Done
897  return ;
898 
899 }

◆ 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 931 of file ProSHADE.cpp.

933 {
934  //================================================ Set the value
935  this->taylorSeriesCap = tayCap;
936 
937  //================================================ Done
938  return ;
939 
940 }

◆ 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 993 of file ProSHADE.cpp.

995 {
996  //================================================ Set the value
997  this->computeTraceSigmaDesc = trSigVal;
998 
999  //================================================ Done
1000  return ;
1001 
1002 }

◆ 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 567 of file ProSHADE.cpp.

569 {
570  //================================================ Set the value
571  this->halfMapKernel = typNoi;
572 
573  //================================================ Done
574  return ;
575 
576 }

◆ 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 463 of file ProSHADE.cpp.

465 {
466  //================================================ Set the value
467  this->verbose = verbosity;
468 
469  //================================================ Done
470  return ;
471 
472 }

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:118
ProSHADE_settings::setOverlayJsonFile
void setOverlayJsonFile(std::string filename)
Sets the filename to which the overlay operations are to be save into.
Definition: ProSHADE.cpp:1367
ProSHADE_settings::maxBandwidth
proshade_unsign maxBandwidth
The bandwidth of spherical harmonics decomposition for the largest sphere.
Definition: ProSHADE_settings.hpp:58
ProSHADE_settings::integOrder
proshade_unsign integOrder
The order required for full Gauss-Legendre integration between the spheres.
Definition: ProSHADE_settings.hpp:68
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:128
ProSHADE_settings::rotTrsJSONFile
std::string rotTrsJSONFile
The filename to which the rotation and translation operations are to be saved into.
Definition: ProSHADE_settings.hpp:139
ProSHADE_settings::setEnLevShellWeight
void setEnLevShellWeight(proshade_double mPower)
Sets the weight of shell position for the energy levels computation.
Definition: ProSHADE.cpp:1101
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:113
ProSHADE_settings::setTraceSigmaComputation
void setTraceSigmaComputation(bool trSigVal)
Sets whether the trace sigma distance descriptor should be computed.
Definition: ProSHADE.cpp:993
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:114
ProSHADE_settings::setSameBoundaries
void setSameBoundaries(bool sameB)
Sets whether same boundaries should be used in the appropriate variable.
Definition: ProSHADE.cpp:749
ProSHADE_settings::maxSymmetryFold
proshade_unsign maxSymmetryFold
The highest symmetry fold to search for.
Definition: ProSHADE_settings.hpp:133
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:69
ProSHADE_settings::setAppliedMaskFilename
void setAppliedMaskFilename(std::string mskFln)
Sets the filename of the mask data that should be applied to the input map.
Definition: ProSHADE.cpp:647
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:93
ProSHADE_settings::setBicubicInterpolationSearch
void setBicubicInterpolationSearch(bool bicubPeaks)
Sets the bicubic interpolation on peaks.
Definition: ProSHADE.cpp:1385
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:1057
ProSHADE_settings::outName
std::string outName
The file name where the output structure(s) should be saved.
Definition: ProSHADE_settings.hpp:108
ProSHADE_settings::setMapInversion
void setMapInversion(bool mInv)
Sets the requested map inversion value in the appropriate variable.
Definition: ProSHADE.cpp:443
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:1035
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:78
ProSHADE_settings::maskMap
bool maskMap
Should the map be masked from noise?
Definition: ProSHADE_settings.hpp:80
ProSHADE_settings::requestedSymmetryFold
proshade_unsign requestedSymmetryFold
The fold of the requested symmetry (only applicable to C and D symmetry types).
Definition: ProSHADE_settings.hpp:131
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:130
ProSHADE_settings::setSphereDistances
void setSphereDistances(proshade_single sphDist)
Sets the requested distance between spheres in the appropriate variable.
Definition: ProSHADE.cpp:890
ProSHADE_settings::setVerbosity
void setVerbosity(proshade_signed verbosity)
Sets the requested verbosity in the appropriate variable.
Definition: ProSHADE.cpp:463
ProSHADE_settings::removeNegativeDensity
bool removeNegativeDensity
Should the negative density be removed from input files?
Definition: ProSHADE_settings.hpp:47
ProSHADE_settings::saveMask
bool saveMask
Should the mask be saved?
Definition: ProSHADE_settings.hpp:84
ProSHADE_settings::requestedResolution
proshade_single requestedResolution
The resolution to which the calculations are to be done.
Definition: ProSHADE_settings.hpp:50
ProSHADE_settings::minSymPeak
proshade_double minSymPeak
Minimum average peak for symmetry axis to be considered as "real".
Definition: ProSHADE_settings.hpp:127
ProSHADE_settings::setAxisComparisonThresholdBehaviour
void setAxisComparisonThresholdBehaviour(bool behav)
Sets the automatic symmetry axis tolerance decreasing.
Definition: ProSHADE.cpp:1188
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:105
ProSHADE_settings::maxSphereDists
proshade_single maxSphereDists
The distance between spheres in spherical mapping for the largest sphere.
Definition: ProSHADE_settings.hpp:65
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:124
ProSHADE_settings::setPDBBFactor
void setPDBBFactor(proshade_double newBF)
Sets the requested B-factor value for PDB files in the appropriate variable.
Definition: ProSHADE.cpp:403
ProSHADE_settings::addStructure
void addStructure(std::string structure)
Adds a structure file name to the appropriate variable.
Definition: ProSHADE.cpp:363
ProSHADE_settings::determineIntegrationOrder
void determineIntegrationOrder(proshade_single maxMapRange)
This function determines the integration order for the between spheres integration.
Definition: ProSHADE.cpp:1571
ProSHADE_settings::setMaskIQR
void setMaskIQR(proshade_single noIQRs)
Sets the requested number of IQRs for masking threshold in the appropriate variable.
Definition: ProSHADE.cpp:504
ProSHADE_settings::setMissingPeakThreshold
void setMissingPeakThreshold(proshade_double mpThres)
Sets the threshold for starting the missing peaks procedure.
Definition: ProSHADE.cpp:1144
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:79
ProSHADE_settings::useCorrelationMasking
bool useCorrelationMasking
Should the blurring masking (false) or the correlation masking (true) be used?
Definition: ProSHADE_settings.hpp:81
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:62
ProSHADE_settings::setMaskBlurFactor
void setMaskBlurFactor(proshade_single blurFac)
Sets the requested map blurring factor in the appropriate variable.
Definition: ProSHADE.cpp:483
ProSHADE_settings::setPeakThreshold
void setPeakThreshold(proshade_double peakThr)
Sets the minimum peak height threshold for axis to be considered possible.
Definition: ProSHADE.cpp:1439
ProSHADE_settings::setPhaseUsage
void setPhaseUsage(bool phaseUsage)
Sets whether the phase information will be used.
Definition: ProSHADE.cpp:1079
ProSHADE_settings::maskFileName
std::string maskFileName
The filename to which mask should be saved.
Definition: ProSHADE_settings.hpp:85
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:142
ProSHADE_settings::setFSCThreshold
void setFSCThreshold(proshade_double fscThr)
Sets the minimum FSC threshold for axis to be considered detected.
Definition: ProSHADE.cpp:1421
ProSHADE_settings::setBandwidth
void setBandwidth(proshade_unsign band)
Sets the requested spherical harmonics bandwidth in the appropriate variable.
Definition: ProSHADE.cpp:870
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:708
ProSHADE_settings::changeMapResolutionTriLinear
bool changeMapResolutionTriLinear
Should maps be re-sampled to obtain the required resolution?
Definition: ProSHADE_settings.hpp:52
ProSHADE_settings::setOutputFilename
void setOutputFilename(std::string oFileName)
Sets the requested output file name in the appropriate variable.
Definition: ProSHADE.cpp:770
ProSHADE_settings::useBiCubicInterpolationOnPeaks
bool useBiCubicInterpolationOnPeaks
This variable switch decides whether best symmetry is detected from peak indices, or whether bicubic ...
Definition: ProSHADE_settings.hpp:132
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:537
ProSHADE_settings::setAxisComparisonThreshold
void setAxisComparisonThreshold(proshade_double axThres)
Sets the threshold for matching symmetry axes.
Definition: ProSHADE.cpp:1165
ProSHADE_settings::peakThresholdMin
proshade_double peakThresholdMin
The threshold for peak height above which axes are considered possible.
Definition: ProSHADE_settings.hpp:135
ProSHADE_settings::forceP1
bool forceP1
Should the P1 spacegroup be forced on the input PDB files?
Definition: ProSHADE_settings.hpp:44
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:310
ProSHADE_settings::task
ProSHADE_Task task
This custom type variable determines which task to perfom (i.e. symmetry detection,...
Definition: ProSHADE_settings.hpp:40
ProSHADE_settings::setMaskFilename
void setMaskFilename(std::string mskFln)
Sets where the mask should be saved.
Definition: ProSHADE.cpp:627
ProSHADE_settings::reBoxMap
bool reBoxMap
This switch decides whether re-boxing is needed.
Definition: ProSHADE_settings.hpp:92
ProSHADE_settings::appliedMaskFileName
std::string appliedMaskFileName
The filename from which mask data will be read from.
Definition: ProSHADE_settings.hpp:86
ProSHADE_settings::normaliseMap
bool normaliseMap
Should the map be normalised to mean 0 sd 1?
Definition: ProSHADE_settings.hpp:72
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:102
ProSHADE_settings::setRequestedFold
void setRequestedFold(proshade_unsign val)
Sets the user requested symmetry fold.
Definition: ProSHADE.cpp:1294
ProSHADE_settings::fscThreshold
proshade_double fscThreshold
The threshold for FSC value under which the axis is considered to be likely noise.
Definition: ProSHADE_settings.hpp:134
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:1508
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:687
ProSHADE_settings::setMaxSymmetryFold
void setMaxSymmetryFold(proshade_unsign maxFold)
Sets the maximum symmetry fold (well, the maximum prime symmetry fold).
Definition: ProSHADE.cpp:1403
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:83
ProSHADE_settings::setNegativeDensity
void setNegativeDensity(bool nDens)
Sets the internal variable deciding whether input files negative density should be removed.
Definition: ProSHADE.cpp:1457
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:75
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:931
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:94
ProSHADE_settings::setProgressiveSphereMapping
void setProgressiveSphereMapping(bool progSphMap)
Sets the requested sphere mapping value settings approach in the appropriate variable.
Definition: ProSHADE.cpp:951
ProSHADE_settings::changeMapResolution
bool changeMapResolution
Should maps be re-sampled to obtain the required resolution?
Definition: ProSHADE_settings.hpp:51
ProSHADE_settings::setIntegrationOrder
void setIntegrationOrder(proshade_unsign intOrd)
Sets the requested order for the Gauss-Legendre integration in the appropriate variable.
Definition: ProSHADE.cpp:910
ProSHADE_settings::setRotationFunctionComputation
void setRotationFunctionComputation(bool rotfVal)
Sets whether the rotation function distance descriptor should be computed.
Definition: ProSHADE.cpp:1014
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:59
ProSHADE_settings::setOverlaySaveFile
void setOverlaySaveFile(std::string filename)
Sets the filename to which the overlay structure is to be save into.
Definition: ProSHADE.cpp:1349
ProSHADE_settings::overlayStructureName
std::string overlayStructureName
The filename to which the rotated and translated moving structure is to be saved.
Definition: ProSHADE_settings.hpp:138
ProSHADE_settings::determineBandwidth
void determineBandwidth(proshade_unsign circumference)
This function determines the bandwidth for the spherical harmonics computation.
Definition: ProSHADE.cpp:1476
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:99
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:117
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:561
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:121
ProSHADE_settings::setFourierWeightsFilename
void setFourierWeightsFilename(std::string fWgFln)
Sets the filename of the mask data that should be applied to the input map.
Definition: ProSHADE.cpp:667
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:82
ProSHADE_settings::determineSphereDistances
void determineSphereDistances(proshade_single maxMapRange)
This function determines the sphere distances for sphere mapping.
Definition: ProSHADE.cpp:1539
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:67
ProSHADE_internal_spheres::autoDetermineBandwidth
proshade_unsign autoDetermineBandwidth(proshade_unsign circumference)
This function determines the bandwidth for the spherical harmonics computation.
Definition: ProSHADE_spheres.cpp:515
ProSHADE_settings::setRequestedSymmetry
void setRequestedSymmetry(std::string val)
Sets the user requested symmetry type.
Definition: ProSHADE.cpp:1274
ProSHADE_settings::setExtraSpace
void setExtraSpace(proshade_single exSpace)
Sets the requested map extra space value in the appropriate variable.
Definition: ProSHADE.cpp:850
ProSHADE_settings::forceBounds
proshade_signed * forceBounds
These will be the boundaries to be forced upon the map.
Definition: ProSHADE_settings.hpp:96
ProSHADE_settings::fourierWeightsFileName
std::string fourierWeightsFileName
The filename from which Fourier weights data will be read from.
Definition: ProSHADE_settings.hpp:89
ProSHADE_settings::detectedSymmetry
std::vector< proshade_double * > detectedSymmetry
The vector of detected symmetry axes.
Definition: ProSHADE_settings.hpp:146
ProSHADE_settings::recommendedSymmetryFold
proshade_unsign recommendedSymmetryFold
The fold of the recommended symmetry C or D type, 0 otherwise.
Definition: ProSHADE_settings.hpp:129
ProSHADE_settings::inputFiles
std::vector< std::string > inputFiles
This vector contains the filenames of all input structure files.
Definition: ProSHADE_settings.hpp:43
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:810
ProSHADE_settings::setMaskSaving
void setMaskSaving(bool savMsk)
Sets whether the mask should be saved.
Definition: ProSHADE.cpp:607
ProSHADE_settings::firstModelOnly
bool firstModelOnly
Shoud only the first PDB model be used, or should all models be used?
Definition: ProSHADE_settings.hpp:46
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:728
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:125
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:95
ProSHADE_settings::setMasking
void setMasking(bool mask)
Sets the requested map masking decision value in the appropriate variable.
Definition: ProSHADE.cpp:524
ProSHADE_settings::pdbBFactorNewVal
proshade_double pdbBFactorNewVal
Change all PDB B-factors to this value (for smooth maps).
Definition: ProSHADE_settings.hpp:55
ProSHADE_settings::removeWaters
bool removeWaters
Should all waters be removed from input PDB files?
Definition: ProSHADE_settings.hpp:45
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:383
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:112
ProSHADE_settings::setNormalisation
void setNormalisation(bool normalise)
Sets the requested map normalisation value in the appropriate variable.
Definition: ProSHADE.cpp:423
ProSHADE_settings::setEnergyLevelsComputation
void setEnergyLevelsComputation(bool enLevDesc)
Sets whether the energy level distance descriptor should be computed.
Definition: ProSHADE.cpp:972
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:111