ProSHADE  0.7.6.0 (JUL 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 setMapReboxing (bool reBx)
 Sets whether re-boxing needs to be done in the appropriate variable. More...
 
void setBoundsSpace (proshade_single boundsExSp)
 Sets the requested number of angstroms for extra space in re-boxing in the appropriate variable. More...
 
void setBoundsThreshold (proshade_signed boundsThres)
 Sets the threshold for number of indices difference acceptable to make index sizes same in the appropriate variable. More...
 
void setSameBoundaries (bool sameB)
 Sets whether same boundaries should be used in the appropriate variable. More...
 
void setOutputFilename (std::string oFileName)
 Sets the requested output file name in the appropriate variable. More...
 
void setMapResolutionChange (bool mrChange)
 Sets the requested map resolution change decision in the appropriate variable. More...
 
void setMapResolutionChangeTriLinear (bool mrChange)
 Sets the requested map resolution change decision using tri-linear interpolation in the appropriate variable. More...
 
void setMapCentering (bool com)
 Sets the requested map centering decision value in the appropriate variable. More...
 
void setExtraSpace (proshade_single exSpace)
 Sets the requested map extra space value in the appropriate variable. More...
 
void setBandwidth (proshade_unsign band)
 Sets the requested spherical harmonics bandwidth in the appropriate variable. More...
 
void setSphereDistances (proshade_single sphDist)
 Sets the requested distance between spheres in the appropriate variable. More...
 
void setIntegrationOrder (proshade_unsign intOrd)
 Sets the requested order for the Gauss-Legendre integration in the appropriate variable. More...
 
void setTaylorSeriesCap (proshade_unsign tayCap)
 Sets the requested Taylor series cap for the Gauss-Legendre integration in the appropriate variable. More...
 
void setProgressiveSphereMapping (bool progSphMap)
 Sets the requested sphere mapping value settings approach in the appropriate variable. More...
 
void setEnergyLevelsComputation (bool enLevDesc)
 Sets whether the energy level distance descriptor should be computed. More...
 
void setTraceSigmaComputation (bool trSigVal)
 Sets whether the trace sigma distance descriptor should be computed. More...
 
void setRotationFunctionComputation (bool rotfVal)
 Sets whether the rotation function distance descriptor should be computed. More...
 
void setPeakNeighboursNumber (proshade_unsign pkS)
 Sets the number of neighbour values that have to be smaller for an index to be considered a peak. More...
 
void setPeakNaiveNoIQR (proshade_double noIQRs)
 Sets the number of IQRs from the median for threshold height a peak needs to be considered a peak. More...
 
void setPhaseUsage (bool phaseUsage)
 Sets whether the phase information will be used. More...
 
void setEnLevShellWeight (proshade_double mPower)
 Sets the weight of shell position for the energy levels computation. More...
 
void setGroupingSmoothingFactor (proshade_double smFact)
 Sets the grouping smoothing factor into the proper variable. More...
 
void setMissingPeakThreshold (proshade_double mpThres)
 Sets the threshold for starting the missing peaks procedure. More...
 
void setAxisComparisonThreshold (proshade_double axThres)
 Sets the threshold for matching symmetry axes. More...
 
void setAxisComparisonThresholdBehaviour (bool behav)
 Sets the automatic symmetry axis tolerance decreasing. More...
 
void setMinimumPeakForAxis (proshade_double minSP)
 Sets the minimum peak height for symmetry axis to be considered. More...
 
void setRecommendedSymmetry (std::string val)
 Sets the ProSHADE detected symmetry type. More...
 
void setRecommendedFold (proshade_unsign val)
 Sets the ProSHADE detected symmetry fold. More...
 
void setRequestedSymmetry (std::string val)
 Sets the user requested symmetry type. More...
 
void setRequestedFold (proshade_unsign val)
 Sets the user requested symmetry fold. More...
 
void setDetectedSymmetry (proshade_double *sym)
 Sets the final detected symmetry axes information. More...
 
void setOverlaySaveFile (std::string filename)
 Sets the filename to which the overlay structure is to be save into. More...
 
void setOverlayJsonFile (std::string filename)
 Sets the filename to which the overlay operations are to be save into. More...
 
void setSymmetryRotFunPeaks (bool rotFunPeaks)
 Sets the symmetry detection algorithm type. More...
 
void setBicubicInterpolationSearch (bool bicubPeaks)
 Sets the bicubic interpolation on peaks. More...
 
void setMaxSymmetryFold (proshade_unsign maxFold)
 Sets the maximum symmetry fold (well, the maximum prime symmetry fold). More...
 
void 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.
 
bool reBoxMap
 This switch decides whether re-boxing is needed.
 
proshade_single boundsExtraSpace
 The number of extra angstroms to be added to all re-boxing bounds just for safety.
 
proshade_signed boundsSimilarityThreshold
 Number of indices which can be added just to make sure same size in indices is achieved.
 
bool useSameBounds
 Switch to say that the same boundaries as used for the first should be used for all input maps.
 
proshade_signed * forceBounds
 These will be the boundaries to be forced upon the map.
 
bool moveToCOM
 Logical value stating whether the structure should be moved to have its Centre Of Mass (COM) in the middle.
 
proshade_single addExtraSpace
 If this value is non-zero, this many angstroms of empty space will be added to the internal map.
 
bool progressiveSphereMapping
 If true, each shell will have its own angular resolution dependent on the actual number of map points which are available to it. If false, all shells will have the same settings.
 
std::string outName
 The file name where the output structure(s) should be saved.
 
bool computeEnergyLevelsDesc
 If true, the energy levels descriptor will be computed, otherwise all its computations will be omitted.
 
proshade_double enLevMatrixPowerWeight
 If RRP matrices shell position is to be weighted by putting the position as an exponent, this variable sets the exponent. Set to 0 for no weighting.
 
bool computeTraceSigmaDesc
 If true, the trace sigma descriptor will be computed, otherwise all its computations will be omitted.
 
bool computeRotationFuncDesc
 If true, the rotation function descriptor will be computed, otherwise all its computations will be omitted.
 
proshade_unsign peakNeighbours
 Number of points in any direction that have to be lower than the considered index in order to consider this index a peak.
 
proshade_double noIQRsFromMedianNaivePeak
 When doing peak searching, how many IQRs from the median the threshold for peak height should be (in terms of median of non-peak values).
 
proshade_double smoothingFactor
 This factor decides how small the group sizes should be - larger factor means more smaller groups.
 
proshade_double symMissPeakThres
 Percentage of peaks that could be missing that would warrant starting the missing peaks search procedure.
 
proshade_double axisErrTolerance
 Allowed error on vector axis in in dot product ( acos ( 1 - axErr ) is the allowed difference in radians ).
 
bool axisErrToleranceDefault
 
proshade_double minSymPeak
 Minimum average peak for symmetry axis to be considered as "real".
 
std::string recommendedSymmetryType
 The symmetry type that ProSHADE finds the best fitting for the structure. Possible values are "" for none, "C" for cyclic, "D" for Dihedral, "T" for Tetrahedral, "O" for Octahedral and "I" for Icosahedral. C and D types also have fold value associated.
 
proshade_unsign recommendedSymmetryFold
 The fold of the recommended symmetry C or D type, 0 otherwise.
 
std::string requestedSymmetryType
 The symmetry type requested by the user. Allowed values are C, D, T, O and I.
 
proshade_unsign requestedSymmetryFold
 The fold of the requested symmetry (only applicable to C and D symmetry types).
 
bool usePeakSearchInRotationFunctionSpace
 This variable switch decides whether symmetry detection will be done using peak search in rotation function or using the angle-axis sperical space.
 
bool useBiCubicInterpolationOnPeaks
 This variable switch decides whether best symmetry is detected from peak indices, or whether bicubic interpolation is done to seatch for better axis between indices.
 
proshade_unsign maxSymmetryFold
 The highest symmetry fold to search for.
 
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 re-boxing
89  this->reBoxMap = false;
90  this->boundsExtraSpace = 3.0;
91  this->boundsSimilarityThreshold = 0;
92  this->useSameBounds = false;
93  this->forceBounds = new proshade_signed [6];
94 
95  //================================================ Settings regarding COM
96  this->moveToCOM = false;
97 
98  //================================================ Settings regarding extra cell space
99  this->addExtraSpace = 10.0;
100 
101  //================================================ Settings regarding shell settings
102  this->progressiveSphereMapping = false;
103 
104  //================================================ Settings regarding output file name
105  this->outName = "reBoxed";
106 
107  //================================================ Settings regarding distances computation
108  this->computeEnergyLevelsDesc = true;
109  this->computeTraceSigmaDesc = true;
110  this->computeRotationFuncDesc = true;
111  this->enLevMatrixPowerWeight = 1.0;
112 
113  //================================================ Settings regarding peak searching
114  this->peakNeighbours = 1;
115  this->noIQRsFromMedianNaivePeak = -999.9;
116 
117  //================================================ Settings regarding 1D grouping
118  this->smoothingFactor = 15.0;
119 
120  //================================================ Settings regarding the symmetry detection
122  this->useBiCubicInterpolationOnPeaks = true;
123  this->maxSymmetryFold = 30;
124  this->fscThreshold = 0.80;
125  this->peakThresholdMin = 0.80;
126  this->symMissPeakThres = 0.3;
127  this->axisErrTolerance = 0.01;
128  this->axisErrToleranceDefault = true;
129  this->minSymPeak = 0.3;
130  this->recommendedSymmetryType = "";
131  this->recommendedSymmetryFold = 0;
132  this->requestedSymmetryType = "";
133  this->requestedSymmetryFold = 0;
134  this->detectedSymmetry.clear ( );
135 
136  //================================================ Settings regarding the structure overlay
137  this->overlayStructureName = "movedStructure";
138  this->rotTrsJSONFile = "movedStructureOperations.json";
139 
140  //================================================ Settings regarding verbosity of the program
141  this->verbose = 1;
142 
143  //================================================ Done
144 
145 }

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

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

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

320 {
321  //================================================ Release boundaries variable
322  delete[] this->forceBounds;
323 
324  //================================================ Release symmetry axes
325  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); } } }
326 
327  //================================================ Done
328 
329 }

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

1169 {
1170  //================================================ Set the value
1171  this->axisErrToleranceDefault = behav;
1172 
1173  //================================================ Done
1174  return ;
1175 
1176 }

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

851 {
852  //================================================ Set the value
853  this->maxBandwidth = band;
854 
855  //================================================ Done
856  return ;
857 
858 }

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

1384 {
1385  //================================================ Set the value
1386  this->useBiCubicInterpolationOnPeaks = bicubPeaks;
1387 
1388  //================================================ Done
1389  return ;
1390 
1391 }

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

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

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

709 {
710  //================================================ Set the value
711  this->boundsSimilarityThreshold = boundsThres;
712 
713  //================================================ Done
714  return ;
715 
716 }

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

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

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

1296 {
1297  //================================================ Allocate memory
1298  proshade_double* hlpAxis = new proshade_double [7];
1299  ProSHADE_internal_misc::checkMemoryAllocation ( hlpAxis, __FILE__, __LINE__, __func__ );
1300 
1301  //================================================ Copy (deep) data
1302  hlpAxis[0] = sym[0];
1303  hlpAxis[1] = sym[1];
1304  hlpAxis[2] = sym[2];
1305  hlpAxis[3] = sym[3];
1306  hlpAxis[4] = sym[4];
1307  hlpAxis[5] = sym[5];
1308  hlpAxis[6] = sym[6];
1309 
1310  //================================================ Save
1312 
1313  //================================================ Release memory
1314  delete[] hlpAxis;
1315 
1316  //================================================ Done
1317  return ;
1318 
1319 }

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

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

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

1082 {
1083  //================================================ Set the value
1084  this->enLevMatrixPowerWeight = mPower;
1085 
1086  //================================================ Done
1087  return ;
1088 
1089 }

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

831 {
832  //================================================ Set the value
833  this->addExtraSpace = exSpace;
834 
835  //================================================ Done
836  return ;
837 
838 }

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

1420 {
1421  //================================================ Set the value
1422  this->fscThreshold = fscThr;
1423 
1424  //================================================ Done
1425  return ;
1426 
1427 }

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

1104 {
1105  //================================================ Set the value
1106  this->smoothingFactor = smFact;
1107 
1108  //================================================ Done
1109  return ;
1110 
1111 }

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

891 {
892  //================================================ Set the value
893  this->integOrder = intOrd;
894 
895  //================================================ Done
896  return ;
897 
898 }

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

811 {
812  //================================================ Set the value
813  this->moveToCOM = com;
814 
815  //================================================ Done
816  return ;
817 
818 }

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

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

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

668 {
669  //================================================ Set the value
670  this->reBoxMap = reBx;
671 
672  //================================================ Done
673  return ;
674 
675 }

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

771 {
772  //================================================ Set the value
773  this->changeMapResolution = mrChange;
774 
775  //================================================ Done
776  return ;
777 
778 }

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

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

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

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

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

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

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

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

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

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

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

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

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

1402 {
1403  //================================================ Set the value
1404  this->maxSymmetryFold = maxFold;
1405 
1406  //================================================ Done
1407  return ;
1408 
1409 }

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

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

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

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

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

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

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

1456 {
1457  //================================================ Set the value
1458  this->removeNegativeDensity = nDens;
1459 
1460  //================================================ Done
1461  return ;
1462 
1463 }

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

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

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

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

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

1348 {
1349  //================================================ Set the value
1350  this->rotTrsJSONFile = filename;
1351 
1352  //================================================ Done
1353  return ;
1354 
1355 }

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

1330 {
1331  //================================================ Set the value
1332  this->overlayStructureName = filename;
1333 
1334  //================================================ Done
1335  return ;
1336 
1337 }

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

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

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

1038 {
1039  //================================================ Set the value
1040  this->noIQRsFromMedianNaivePeak = noIQRs;
1041 
1042  //================================================ Done
1043  return ;
1044 
1045 }

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

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

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

1438 {
1439  //================================================ Set the value
1440  this->peakThresholdMin = peakThr;
1441 
1442  //================================================ Done
1443  return ;
1444 
1445 }

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

1060 {
1061  //================================================ Set the value
1062  this->usePhase = phaseUsage;
1063 
1064  //================================================ Done
1065  return ;
1066 
1067 }

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

932 {
933  //================================================ Set the value
934  this->progressiveSphereMapping = progSphMap;
935 
936  //================================================ Done
937  return ;
938 
939 }

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

1235 {
1236  //================================================ Set the value
1237  this->recommendedSymmetryFold = val;
1238 
1239  //================================================ Done
1240  return ;
1241 
1242 }

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

1212 {
1213  //================================================ Set the value
1214  this->recommendedSymmetryType = val;
1215 
1216  //================================================ Done
1217  return ;
1218 
1219 }

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

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

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

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

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

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

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

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

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

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

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

871 {
872  //================================================ Set the value
873  this->maxSphereDists = sphDist;
874 
875  //================================================ Done
876  return ;
877 
878 }

◆ setSymmetryRotFunPeaks()

void ProSHADE_settings::setSymmetryRotFunPeaks ( bool  rotFunPeaks)

Sets the symmetry detection algorithm type.

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

Definition at line 1364 of file ProSHADE.cpp.

1366 {
1367  //================================================ Set the value
1368  this->usePeakSearchInRotationFunctionSpace = rotFunPeaks;
1369 
1370  //================================================ Done
1371  return ;
1372 
1373 }

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

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

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

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

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

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

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

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

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:115
ProSHADE_settings::setOverlayJsonFile
void setOverlayJsonFile(std::string filename)
Sets the filename to which the overlay operations are to be save into.
Definition: ProSHADE.cpp:1346
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:125
ProSHADE_settings::rotTrsJSONFile
std::string rotTrsJSONFile
The filename to which the rotation and translation operations are to be saved into.
Definition: ProSHADE_settings.hpp:137
ProSHADE_settings::setEnLevShellWeight
void setEnLevShellWeight(proshade_double mPower)
Sets the weight of shell position for the energy levels computation.
Definition: ProSHADE.cpp:1080
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:110
ProSHADE_settings::setTraceSigmaComputation
void setTraceSigmaComputation(bool trSigVal)
Sets whether the trace sigma distance descriptor should be computed.
Definition: ProSHADE.cpp:972
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:111
ProSHADE_settings::setSameBoundaries
void setSameBoundaries(bool sameB)
Sets whether same boundaries should be used in the appropriate variable.
Definition: ProSHADE.cpp:728
ProSHADE_settings::maxSymmetryFold
proshade_unsign maxSymmetryFold
The highest symmetry fold to search for.
Definition: ProSHADE_settings.hpp:131
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:646
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:90
ProSHADE_settings::setBicubicInterpolationSearch
void setBicubicInterpolationSearch(bool bicubPeaks)
Sets the bicubic interpolation on peaks.
Definition: ProSHADE.cpp:1382
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:1036
ProSHADE_settings::outName
std::string outName
The file name where the output structure(s) should be saved.
Definition: ProSHADE_settings.hpp:105
ProSHADE_settings::setMapInversion
void setMapInversion(bool mInv)
Sets the requested map inversion value in the appropriate variable.
Definition: ProSHADE.cpp:442
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:1014
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:128
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:127
ProSHADE_settings::setSphereDistances
void setSphereDistances(proshade_single sphDist)
Sets the requested distance between spheres in the appropriate variable.
Definition: ProSHADE.cpp:869
ProSHADE_settings::setVerbosity
void setVerbosity(proshade_signed verbosity)
Sets the requested verbosity in the appropriate variable.
Definition: ProSHADE.cpp:462
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:124
ProSHADE_settings::setAxisComparisonThresholdBehaviour
void setAxisComparisonThresholdBehaviour(bool behav)
Sets the automatic symmetry axis tolerance decreasing.
Definition: ProSHADE.cpp:1167
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:102
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:121
ProSHADE_settings::setPDBBFactor
void setPDBBFactor(proshade_double newBF)
Sets the requested B-factor value for PDB files in the appropriate variable.
Definition: ProSHADE.cpp:402
ProSHADE_settings::addStructure
void addStructure(std::string structure)
Adds a structure file name to the appropriate variable.
Definition: ProSHADE.cpp:362
ProSHADE_settings::determineIntegrationOrder
void determineIntegrationOrder(proshade_single maxMapRange)
This function determines the integration order for the between spheres integration.
Definition: ProSHADE.cpp:1568
ProSHADE_settings::setMaskIQR
void setMaskIQR(proshade_single noIQRs)
Sets the requested number of IQRs for masking threshold in the appropriate variable.
Definition: ProSHADE.cpp:503
ProSHADE_settings::setMissingPeakThreshold
void setMissingPeakThreshold(proshade_double mpThres)
Sets the threshold for starting the missing peaks procedure.
Definition: ProSHADE.cpp:1123
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:482
ProSHADE_settings::setPeakThreshold
void setPeakThreshold(proshade_double peakThr)
Sets the minimum peak height threshold for axis to be considered possible.
Definition: ProSHADE.cpp:1436
ProSHADE_settings::setPhaseUsage
void setPhaseUsage(bool phaseUsage)
Sets whether the phase information will be used.
Definition: ProSHADE.cpp:1058
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:140
ProSHADE_settings::setFSCThreshold
void setFSCThreshold(proshade_double fscThr)
Sets the minimum FSC threshold for axis to be considered detected.
Definition: ProSHADE.cpp:1418
ProSHADE_settings::setBandwidth
void setBandwidth(proshade_unsign band)
Sets the requested spherical harmonics bandwidth in the appropriate variable.
Definition: ProSHADE.cpp:849
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:687
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:749
ProSHADE_settings::useBiCubicInterpolationOnPeaks
bool useBiCubicInterpolationOnPeaks
This variable switch decides whether best symmetry is detected from peak indices, or whether bicubic ...
Definition: ProSHADE_settings.hpp:130
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:1144
ProSHADE_settings::peakThresholdMin
proshade_double peakThresholdMin
The threshold for peak height above which axes are considered possible.
Definition: ProSHADE_settings.hpp:133
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:626
ProSHADE_settings::reBoxMap
bool reBoxMap
This switch decides whether re-boxing is needed.
Definition: ProSHADE_settings.hpp:89
ProSHADE_settings::setSymmetryRotFunPeaks
void setSymmetryRotFunPeaks(bool rotFunPeaks)
Sets the symmetry detection algorithm type.
Definition: ProSHADE.cpp:1364
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:99
ProSHADE_settings::setRequestedFold
void setRequestedFold(proshade_unsign val)
Sets the user requested symmetry fold.
Definition: ProSHADE.cpp:1273
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:132
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:1505
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:666
ProSHADE_settings::setMaxSymmetryFold
void setMaxSymmetryFold(proshade_unsign maxFold)
Sets the maximum symmetry fold (well, the maximum prime symmetry fold).
Definition: ProSHADE.cpp:1400
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:1454
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:910
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:91
ProSHADE_settings::setProgressiveSphereMapping
void setProgressiveSphereMapping(bool progSphMap)
Sets the requested sphere mapping value settings approach in the appropriate variable.
Definition: ProSHADE.cpp:930
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:889
ProSHADE_settings::setRotationFunctionComputation
void setRotationFunctionComputation(bool rotfVal)
Sets whether the rotation function distance descriptor should be computed.
Definition: ProSHADE.cpp:993
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:1328
ProSHADE_settings::overlayStructureName
std::string overlayStructureName
The filename to which the rotated and translated moving structure is to be saved.
Definition: ProSHADE_settings.hpp:136
ProSHADE_settings::determineBandwidth
void determineBandwidth(proshade_unsign circumference)
This function determines the bandwidth for the spherical harmonics computation.
Definition: ProSHADE.cpp:1473
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:96
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:114
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:118
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:1536
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:1253
ProSHADE_settings::setExtraSpace
void setExtraSpace(proshade_single exSpace)
Sets the requested map extra space value in the appropriate variable.
Definition: ProSHADE.cpp:829
ProSHADE_settings::forceBounds
proshade_signed * forceBounds
These will be the boundaries to be forced upon the map.
Definition: ProSHADE_settings.hpp:93
ProSHADE_settings::detectedSymmetry
std::vector< proshade_double * > detectedSymmetry
The vector of detected symmetry axes.
Definition: ProSHADE_settings.hpp:144
ProSHADE_settings::recommendedSymmetryFold
proshade_unsign recommendedSymmetryFold
The fold of the recommended symmetry C or D type, 0 otherwise.
Definition: ProSHADE_settings.hpp:126
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:789
ProSHADE_settings::setMaskSaving
void setMaskSaving(bool savMsk)
Sets whether the mask should be saved.
Definition: ProSHADE.cpp:606
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:707
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:122
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:92
ProSHADE_settings::setMasking
void setMasking(bool mask)
Sets the requested map masking decision value in the appropriate variable.
Definition: ProSHADE.cpp:523
ProSHADE_settings::pdbBFactorNewVal
proshade_double pdbBFactorNewVal
Change all PDB B-factors to this value (for smooth maps).
Definition: ProSHADE_settings.hpp:55
ProSHADE_settings::usePeakSearchInRotationFunctionSpace
bool usePeakSearchInRotationFunctionSpace
This variable switch decides whether symmetry detection will be done using peak search in rotation fu...
Definition: ProSHADE_settings.hpp:129
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:382
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:109
ProSHADE_settings::setNormalisation
void setNormalisation(bool normalise)
Sets the requested map normalisation value in the appropriate variable.
Definition: ProSHADE.cpp:422
ProSHADE_settings::setEnergyLevelsComputation
void setEnergyLevelsComputation(bool enLevDesc)
Sets whether the energy level distance descriptor should be computed.
Definition: ProSHADE.cpp:951
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:108