ProSHADE  0.7.6.0 (JUL 2021)
Protein Shape Detection
ProSHADE_tasks.cpp
Go to the documentation of this file.
1 
22 //==================================================== ProSHADE
23 #include "ProSHADE_tasks.hpp"
24 
35 void ProSHADE_internal_tasks::MapManipulationTask ( ProSHADE_settings* settings, std::vector < proshade_signed* >* originalBounds, std::vector < proshade_signed* >* reboxedBounds, std::vector < proshade_double* >* manipulatedMaps )
36 {
37  //================================================ Check the settings are complete and meaningful
38  checkMapManipulationSettings ( settings );
39 
40  //================================================ For all inputted structures
41  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( settings->inputFiles.size() ); iter++ )
42  {
43  //============================================ Create a data object
45 
46  //============================================ Read in the file
47  strToRebox->readInStructure ( settings->inputFiles.at(iter), iter, settings );
48 
49  //============================================ Save the original boundaries
50  ProSHADE_internal_misc::deepCopyBoundsSigPtrVector ( originalBounds, strToRebox->getXFromPtr(), strToRebox->getXToPtr(), strToRebox->getYFromPtr(), strToRebox->getYToPtr(), strToRebox->getZFromPtr(), strToRebox->getZToPtr() );
51 
52  //============================================ Internal data processing (COM, norm, mask, extra space)
53  strToRebox->processInternalMap ( settings );
54 
55  //============================================ Create new structure for re-boxing
57 
58  //============================================ Re-box map, if need be
59  if ( settings->reBoxMap )
60  {
61  //======================================== Find non-zero bounds
62  proshade_signed* nonZeroBounds = new proshade_signed[6];
63  strToRebox->getReBoxBoundaries ( settings, nonZeroBounds );
64 
65  //============================================ Create new structure from the bounds
66  strToRebox->createNewMapFromBounds ( settings, reBoxStr, nonZeroBounds );
67 
68  //======================================== Release memory
69  delete[] nonZeroBounds;
70  }
71 
72  //============================================ Save the modified structure
73  std::stringstream ss;
74  ss << settings->outName << "_" << iter << ".map";
75  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 1, "Saving the re-boxed map into " + ss.str() );
76  if ( settings->reBoxMap ) { reBoxStr->writeMap ( ss.str() ); }
77  else { strToRebox->writeMap ( ss.str() ); }
78  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, "Structure saved." );
79 
80  //============================================ Save the re-boxed boundaries
81  ProSHADE_internal_misc::deepCopyBoundsSigPtrVector ( reboxedBounds, reBoxStr->getXFromPtr(), reBoxStr->getXToPtr(), reBoxStr->getYFromPtr(),
82  reBoxStr->getYToPtr(), reBoxStr->getZFromPtr(), reBoxStr->getZToPtr() );
83 
84  //============================================ Save the map
85  proshade_double* mapCopy = nullptr;
86  reBoxStr->deepCopyMap ( mapCopy, settings->verbose );
87  ProSHADE_internal_misc::addToDblPtrVector ( manipulatedMaps, mapCopy );
88 
89  //============================================ Release memory
90  delete strToRebox;
91  delete reBoxStr;
92  }
93 
94  //================================================ Done
95  return ;
96 
97 }
98 
107 {
108  //================================================ Is there a single file for processing?
109  if ( settings->inputFiles.size () == 0 )
110  {
111  throw ProSHADE_exception ( "There is no input structure for map manipulation.", "EB00002", __FILE__, __LINE__, __func__, "The ProSHADE_settings object does not contain any\n : structure that could be manipulated. Please supply exactly\n : one structure using the addStructure() function." );
112  }
113 
114  //================================================ Is the file type MAP? Warning if not
115  if ( ProSHADE_internal_io::isFilePDB ( settings->inputFiles.at(0) ) )
116  {
117  ProSHADE_internal_messages::printWarningMessage ( settings->verbose, "!!! ProSHADE WARNING !!! The input file is not of the MAP (MRC) format. Will output re-boxed map, but beware that this is simple PDB->MAP conversion and REFMAC5 should be used to compute more appropriate maps.", "WB00004" );
118 
119  //============================================ No resolution for PDB? Problem...
120  if ( settings->requestedResolution == 0.0f )
121  {
122  throw ProSHADE_exception ( "No resolution given for PDB file re-boxing.", "EB00011", __FILE__, __LINE__, __func__, "The ProSHADE_settings object does not contain any\n : resolution value. However, resolution is required when\n : re-boxing structures read from PDB files. Please supply\n : the resolution value using the setResolution() function." );
123  }
124  }
125 
126  //================================================ Is there output file name?
127  if ( settings->outName == "" )
128  {
129  throw ProSHADE_exception ( "No output file name.", "EB00016", __FILE__, __LINE__, __func__, "There is no output file name set in the settings object.\n : Please supply the file name to where the re-boxed map\n : should be saved using the setOutputFilename() function." );
130  }
131 
132  //================================================ Done
133  return ;
134 
135 }
136 
147 void ProSHADE_internal_tasks::DistancesComputationTask ( ProSHADE_settings* settings, std::vector< proshade_double >* enLevs, std::vector< proshade_double >* trSigm, std::vector< proshade_double >* rotFun )
148 {
149  //================================================ Check the settings are complete and meaningful
150  checkDistancesSettings ( settings );
151 
152  //================================================ Create a data object
154 
155  //================================================ Read in the structure all others will be compared to
156  compareAgainst->readInStructure ( settings->inputFiles.at(0), 0, settings );
157 
158  //================================================ Internal data processing (COM, norm, mask, extra space)
159  compareAgainst->processInternalMap ( settings );
160 
161  //================================================ Map to sphere
162  compareAgainst->mapToSpheres ( settings );
163 
164  //================================================ Get spherical harmonics
165  compareAgainst->computeSphericalHarmonics ( settings );
166 
167  //================================================ Now, for each other structure
168  for ( proshade_unsign iter = 1; iter < static_cast<proshade_unsign> ( settings->inputFiles.size() ); iter++ )
169  {
170  //============================================ Create a data object
172 
173  //============================================ Read in the compared structure
174  compareChanging->readInStructure ( settings->inputFiles.at(iter), iter, settings );
175 
176  //============================================ Internal data processing (COM, norm, mask, extra space)
177  compareChanging->processInternalMap ( settings );
178 
179  //============================================ Map to sphere
180  compareChanging->mapToSpheres ( settings );
181 
182  //============================================ Get spherical harmonics
183  compareChanging->computeSphericalHarmonics ( settings );
184 
185  //============================================ Get distances
186  proshade_double enLevDist = 0.0;
187  if ( settings->computeEnergyLevelsDesc ) { enLevDist = ProSHADE_internal_distances::computeEnergyLevelsDescriptor ( compareAgainst, compareChanging, settings ); }
188  else { ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 1, "Energy levels distance computation not required." ); }
189 
190  proshade_double trSigmDist = 0.0;
191  if ( settings->computeTraceSigmaDesc ) { trSigmDist = ProSHADE_internal_distances::computeTraceSigmaDescriptor ( compareAgainst, compareChanging, settings ); }
192  else { ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 1, "Trace sigma distance computation not required." ); }
193 
194  proshade_double rotFunDist = 0.0;
195  if ( settings->computeRotationFuncDesc ) { rotFunDist = ProSHADE_internal_distances::computeRotationunctionDescriptor ( compareAgainst, compareChanging, settings ); }
196  else { ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 1, "Rotation function distance computation not required." ); }
197 
198  //============================================ Save results to the run object
199  ProSHADE_internal_misc::addToDoubleVector ( enLevs, enLevDist );
200  ProSHADE_internal_misc::addToDoubleVector ( trSigm, trSigmDist );
201  ProSHADE_internal_misc::addToDoubleVector ( rotFun, rotFunDist );
202 
203  //============================================ Report results
204  ReportDistancesResults ( settings, settings->inputFiles.at(0), settings->inputFiles.at(iter), enLevDist, trSigmDist, rotFunDist );
205 
206  //============================================ Release the memory
207  delete compareChanging;
208  }
209 
210 
211  //================================================ Release memory
212  delete compareAgainst;
213 
214  //================================================ Done
215  return ;
216 
217 }
218 
228 void ProSHADE_internal_tasks::ReportDistancesResults ( ProSHADE_settings* settings, std::string str1, std::string str2, proshade_double enLevDist, proshade_double trSigmDist, proshade_double rotFunDist )
229 {
230  std::stringstream hlpSS;
231  hlpSS << "Distances between " << str1 << " and " << str2;
232  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 0, hlpSS.str() );
233 
234  std::stringstream hlpSSE;
235  hlpSSE << "Energy levels distance : " << enLevDist;
236  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 0, hlpSSE.str() );
237 
238  std::stringstream hlpSSS;
239  hlpSSS << "Trace sigma distance : " << trSigmDist;
240  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 0, hlpSSS.str() );
241 
242  std::stringstream hlpSSR;
243  hlpSSR << "Rotation function distance: " << rotFunDist;
244  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 0, hlpSSR.str() );
245 
246  //================================================ Done
247  return ;
248 
249 }
250 
259 {
260  //================================================ Are there at least two structures?
261  if ( settings->inputFiles.size () < 2 )
262  {
263  throw ProSHADE_exception ( "There are not enough structures for distance computation.", "ED00012", __FILE__, __LINE__, __func__, "There needs to be at least two structures between which\n : distances are computed. The ProSHADE_settings object\n : contains less than two structures and therefore cannot\n : proceed. Please supply at least two structures by\n : repeatedly using the addStructure() function." );
264  }
265 
266  //================================================ Is there resolution value set?
267  const FloatingPoint< proshade_single > lhs ( settings->requestedResolution ), rhs ( -1.0f );
268  if ( lhs.AlmostEquals ( rhs ) )
269  {
270  throw ProSHADE_exception ( "Resolution value not set.", "ED00013", __FILE__, __LINE__, __func__, "The resolution value was not set. Please set the\n : resolution value for the distance computation by using\n : the setResolution() function." );
271  }
272 
273  //================================================ Done
274  return ;
275 
276 }
277 
287 void ProSHADE_internal_tasks::SymmetryDetectionTask ( ProSHADE_settings* settings, std::vector< proshade_double* >* axes, std::vector < std::vector< proshade_double > >* allCs, std::vector< proshade_double >* mapCOMShift )
288 {
289  //================================================ Check the settings are complete and meaningful
290  checkSymmetrySettings ( settings );
291 
292  //================================================ Now, for each other structure
293  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( settings->inputFiles.size() ); iter++ )
294  {
295  //============================================ Create a data object
297 
298  //============================================ Read in the compared structure
299  symmetryStructure->readInStructure ( settings->inputFiles.at(iter), iter, settings );
300 
301  //============================================ Internal data processing (COM, norm, mask, extra space)
302  symmetryStructure->processInternalMap ( settings );
303 
304  //============================================ Map to sphere
305  symmetryStructure->mapToSpheres ( settings );
306 
307  //============================================ Get spherical harmonics
308  symmetryStructure->computeSphericalHarmonics ( settings );
309 
310  //============================================ Compute auto-rotation map
311  symmetryStructure->computeRotationFunction ( settings );
312 
313  if ( settings->usePeakSearchInRotationFunctionSpace )
314  {
315  //======================================== Detect point groups in the angle-axis space
316  symmetryStructure->detectSymmetryFromAngleAxisSpace ( settings, axes, allCs );
317  }
318  else
319  {
320  //======================================== Detect symmetry using the peak detection in rotation function space
321  symmetryStructure->detectSymmetryInStructure ( settings, axes, allCs );
322  }
323 
324  //============================================ Report results
325  symmetryStructure->reportSymmetryResults ( settings );
326 
327  //============================================ Save internal map shift to run object,
328  ProSHADE_internal_misc::addToDoubleVector ( mapCOMShift, symmetryStructure->mapCOMProcessChangeX );
329  ProSHADE_internal_misc::addToDoubleVector ( mapCOMShift, symmetryStructure->mapCOMProcessChangeY );
330  ProSHADE_internal_misc::addToDoubleVector ( mapCOMShift, symmetryStructure->mapCOMProcessChangeZ );
331 
332  //============================================ Release memory
333  delete symmetryStructure;
334  }
335 
336  //================================================ Done
337  return ;
338 
339 }
340 
349 {
350  //================================================ Are the any structures?
351  if ( settings->inputFiles.size () < 1 )
352  {
353  throw ProSHADE_exception ( "There are not enough structures for symmetry detection.", "ES00028", __FILE__, __LINE__, __func__, "There needs to be at least one structure for which\n : symmetry is to be detected. Please supply at least one\n : structure by using the addStructure() function." );
354  }
355 
356  //================================================ Is the axis tolerance set properly?
357  if ( settings->axisErrTolerance < 0.0 )
358  {
359  throw ProSHADE_exception ( "Symmetry axis detection tolerance set to negative value.", "ES00053", __FILE__, __LINE__, __func__, "The symmetry axis detection tolerance was manually set to\n : negative value. This makes no sense, please supply\n : value >= 0.0." );
360  }
361 
362  //================================================ Done
363  return ;
364 
365 }
366 
377 void ProSHADE_internal_tasks::MapOverlayTask ( ProSHADE_settings* settings, std::vector < proshade_double >* rotationCentre, std::vector < proshade_double >* eulerAngles, std::vector < proshade_double >* finalTranslation )
378 {
379  //================================================ Check the settings are complete and meaningful
380  checkOverlaySettings ( settings );
381 
382  //================================================ Initialise variables
383  proshade_double eulA, eulB, eulG, trsX, trsY, trsZ;
384 
385  //================================================ Create the data objects initially (this time without phase)
388 
389  //================================================ First, run without phase and find best rotation angles
390  settings->usePhase = false;
391  ProSHADE_internal_overlay::getOptimalRotation ( settings, staticStructure, movingStructure, &eulA, &eulB, &eulG );
392 
393  //================================================ Release memory
394  delete staticStructure;
395  delete movingStructure;
396 
397  //================================================ Create the data objects again (this time with phase)
398  staticStructure = new ProSHADE_internal_data::ProSHADE_data ( );
399  movingStructure = new ProSHADE_internal_data::ProSHADE_data ( );
400 
401  //================================================ Now, run with phase and find optimal translation
402  settings->usePhase = true;
403  settings->changeMapResolution = true;
404  ProSHADE_internal_overlay::getOptimalTranslation ( settings, staticStructure, movingStructure, &trsX, &trsY, &trsZ, eulA, eulB, eulG );
405 
406  //================================================ Compute the proper translations using the translation function output
407  ProSHADE_internal_misc::addToDoubleVector ( rotationCentre, movingStructure->originalPdbRotCenX );
408  ProSHADE_internal_misc::addToDoubleVector ( rotationCentre, movingStructure->originalPdbRotCenY );
409  ProSHADE_internal_misc::addToDoubleVector ( rotationCentre, movingStructure->originalPdbRotCenZ );
410  ProSHADE_internal_misc::addToDoubleVector ( finalTranslation, movingStructure->originalPdbTransX );
411  ProSHADE_internal_misc::addToDoubleVector ( finalTranslation, movingStructure->originalPdbTransY );
412  ProSHADE_internal_misc::addToDoubleVector ( finalTranslation, movingStructure->originalPdbTransZ );
413 
414  //================================================ Write out everything
415  movingStructure->writeOutOverlayFiles ( settings, eulA, eulB, eulG, rotationCentre, finalTranslation );
416 
417  //================================================ Save the rotation and rest of translations
418  ProSHADE_internal_misc::addToDoubleVector ( eulerAngles, eulA );
419  ProSHADE_internal_misc::addToDoubleVector ( eulerAngles, eulB );
420  ProSHADE_internal_misc::addToDoubleVector ( eulerAngles, eulG );
421 
422  //================================================ Report results to user
423  movingStructure->reportOverlayResults ( settings, rotationCentre, eulerAngles, finalTranslation );
424 
425  //================================================ Release memory
426  delete staticStructure;
427  delete movingStructure;
428 
429  //================================================ Done
430  return ;
431 
432 }
433 
442 {
443  //================================================ Are the any structures?
444  if ( settings->inputFiles.size () != 2 )
445  {
446  throw ProSHADE_exception ( "There are not enough structures for map overlay\n : computation.", "EO00033", __FILE__, __LINE__, __func__, "There needs to be exactly two structures for map overlay\n : mode to work; the first structure is the static and the\n : second is the moving structure." );
447  }
448 
449  //================================================ Done
450  return ;
451 
452 }
ProSHADE_internal_tasks::MapOverlayTask
void MapOverlayTask(ProSHADE_settings *settings, std::vector< proshade_double > *rotationCentre, std::vector< proshade_double > *eulerAngles, std::vector< proshade_double > *finalTranslation)
The symmetry detection task driver function.
Definition: ProSHADE_tasks.cpp:377
ProSHADE_internal_data::ProSHADE_data::computeSphericalHarmonics
void computeSphericalHarmonics(ProSHADE_settings *settings)
This function computes the spherical harmonics decomposition for the whole structure.
Definition: ProSHADE_data.cpp:1756
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::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_internal_misc::addToDblPtrVector
void addToDblPtrVector(std::vector< proshade_double * > *vecToAddTo, proshade_double *elementToAdd)
Adds the element to the vector.
Definition: ProSHADE_misc.cpp:143
ProSHADE_internal_data::ProSHADE_data::createNewMapFromBounds
void createNewMapFromBounds(ProSHADE_settings *settings, ProSHADE_data *&newStr, proshade_signed *newBounds)
This function creates a new structure from the calling structure and new bounds values.
Definition: ProSHADE_data.cpp:1300
ProSHADE_internal_io::isFilePDB
bool isFilePDB(std::string fName)
Function determining if the input data type is PDB.
Definition: ProSHADE_io.cpp:32
ProSHADE_internal_data::ProSHADE_data::computeRotationFunction
void computeRotationFunction(ProSHADE_settings *settings)
This function computes the self-rotation function for this structure.
Definition: ProSHADE_symmetry.cpp:41
ProSHADE_internal_data::ProSHADE_data::getReBoxBoundaries
void getReBoxBoundaries(ProSHADE_settings *settings, proshade_signed *&ret)
This function finds the boundaries enclosing positive map values and adds some extra space.
Definition: ProSHADE_data.cpp:1242
ProSHADE_internal_tasks::checkDistancesSettings
void checkDistancesSettings(ProSHADE_settings *settings)
The distances computation settings checks.
Definition: ProSHADE_tasks.cpp:258
ProSHADE_exception
This class is the representation of ProSHADE exception.
Definition: ProSHADE_exceptions.hpp:37
ProSHADE_settings::outName
std::string outName
The file name where the output structure(s) should be saved.
Definition: ProSHADE_settings.hpp:105
ProSHADE_internal_data::ProSHADE_data
This class contains all inputed and derived data for a single structure.
Definition: ProSHADE_data.hpp:49
ProSHADE_settings::requestedResolution
proshade_single requestedResolution
The resolution to which the calculations are to be done.
Definition: ProSHADE_settings.hpp:50
ProSHADE_internal_data::ProSHADE_data::mapCOMProcessChangeX
proshade_double mapCOMProcessChangeX
The change in X axis between the creation of the structure (originalMapXCom) and just before rotation...
Definition: ProSHADE_data.hpp:97
ProSHADE_internal_data::ProSHADE_data::originalPdbRotCenY
proshade_double originalPdbRotCenY
The centre of rotation as it relates to the original PDB positions (and not the ProSHADE internal map...
Definition: ProSHADE_data.hpp:103
ProSHADE_internal_tasks::checkOverlaySettings
void checkOverlaySettings(ProSHADE_settings *settings)
The map overlay computation settings checks.
Definition: ProSHADE_tasks.cpp:441
ProSHADE_internal_data::ProSHADE_data::mapCOMProcessChangeZ
proshade_double mapCOMProcessChangeZ
The change in Z axis between the creation of the structure (originalMapZCom) and just before rotation...
Definition: ProSHADE_data.hpp:99
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_internal_overlay::getOptimalRotation
void getOptimalRotation(ProSHADE_settings *settings, ProSHADE_internal_data::ProSHADE_data *staticStructure, ProSHADE_internal_data::ProSHADE_data *movingStructure, proshade_double *eulA, proshade_double *eulB, proshade_double *eulG)
This function finds the optimal rotation between two structures as described by the settings object.
Definition: ProSHADE_overlay.cpp:154
ProSHADE_internal_messages::printWarningMessage
void printWarningMessage(proshade_signed verbose, std::string message, std::string warnCode)
General stderr message printing (used for warnings).
Definition: ProSHADE_messages.cpp:101
ProSHADE_internal_distances::computeRotationunctionDescriptor
proshade_double computeRotationunctionDescriptor(ProSHADE_internal_data::ProSHADE_data *obj1, ProSHADE_internal_data::ProSHADE_data *obj2, ProSHADE_settings *settings)
This function computes the rotation function descriptor value between two objects.
Definition: ProSHADE_distances.cpp:906
ProSHADE_internal_data::ProSHADE_data::getXFromPtr
proshade_signed * getXFromPtr(void)
This function allows access to the map start along the X axis.
Definition: ProSHADE_data.cpp:3969
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_internal_data::ProSHADE_data::deepCopyMap
void deepCopyMap(proshade_double *&saveTo, proshade_signed verbose)
This function copies the internal map into the supplied pointer, which it also allocates.
Definition: ProSHADE_data.cpp:3466
ProSHADE_internal_misc::addToDoubleVector
void addToDoubleVector(std::vector< proshade_double > *vecToAddTo, proshade_double elementToAdd)
Adds the element to the vector.
Definition: ProSHADE_misc.cpp:77
ProSHADE_internal_data::ProSHADE_data::originalPdbTransX
proshade_double originalPdbTransX
The optimal translation vector as it relates to the original PDB positions (and not the ProSHADE inte...
Definition: ProSHADE_data.hpp:105
ProSHADE_internal_data::ProSHADE_data::mapToSpheres
void mapToSpheres(ProSHADE_settings *settings)
This function converts the internal map onto a set of concentric spheres.
Definition: ProSHADE_data.cpp:1712
ProSHADE_internal_data::ProSHADE_data::originalPdbTransY
proshade_double originalPdbTransY
The optimal translation vector as it relates to the original PDB positions (and not the ProSHADE inte...
Definition: ProSHADE_data.hpp:106
ProSHADE_internal_tasks::MapManipulationTask
void MapManipulationTask(ProSHADE_settings *settings, std::vector< proshade_signed * > *originalBounds, std::vector< proshade_signed * > *reboxedBounds, std::vector< proshade_double * > *manipulatedMaps)
The re-boxing task driver function.
Definition: ProSHADE_tasks.cpp:35
ProSHADE_internal_data::ProSHADE_data::reportOverlayResults
void reportOverlayResults(ProSHADE_settings *settings, std::vector< proshade_double > *rotationCentre, std::vector< proshade_double > *eulerAngles, std::vector< proshade_double > *finalTranslation)
This function reports the results of the overlay mode.
Definition: ProSHADE_data.cpp:4522
ProSHADE_internal_data::ProSHADE_data::writeMap
void writeMap(std::string fName, std::string title="Created by ProSHADE and written by GEMMI", int mode=2)
Function for writing out the internal structure representation in MRC MAP format.
Definition: ProSHADE_data.cpp:948
ProSHADE_settings::reBoxMap
bool reBoxMap
This switch decides whether re-boxing is needed.
Definition: ProSHADE_settings.hpp:89
ProSHADE_internal_overlay::getOptimalTranslation
void getOptimalTranslation(ProSHADE_settings *settings, ProSHADE_internal_data::ProSHADE_data *staticStructure, ProSHADE_internal_data::ProSHADE_data *movingStructure, proshade_double *trsX, proshade_double *trsY, proshade_double *trsZ, proshade_double eulA, proshade_double eulB, proshade_double eulG)
This function finds the optimal translation between two structures as described by the settings objec...
Definition: ProSHADE_overlay.cpp:203
ProSHADE_settings
This class stores all the settings and is passed to the executive classes instead of a multitude of p...
Definition: ProSHADE_settings.hpp:37
ProSHADE_internal_data::ProSHADE_data::originalPdbTransZ
proshade_double originalPdbTransZ
The optimal translation vector as it relates to the original PDB positions (and not the ProSHADE inte...
Definition: ProSHADE_data.hpp:107
ProSHADE_internal_data::ProSHADE_data::writeOutOverlayFiles
void writeOutOverlayFiles(ProSHADE_settings *settings, proshade_double eulA, proshade_double eulB, proshade_double eulG, std::vector< proshade_double > *rotCentre, std::vector< proshade_double > *ultimateTranslation)
This function writes out the rotated map, co-ordinates and transformation JSON file.
Definition: ProSHADE_data.cpp:4488
ProSHADE_internal_tasks::DistancesComputationTask
void DistancesComputationTask(ProSHADE_settings *settings, std::vector< proshade_double > *enLevs, std::vector< proshade_double > *trSigm, std::vector< proshade_double > *rotFun)
The distances computation task driver function.
Definition: ProSHADE_tasks.cpp:147
ProSHADE_internal_data::ProSHADE_data::processInternalMap
void processInternalMap(ProSHADE_settings *settings)
This function simply clusters several other functions which should be called together.
Definition: ProSHADE_data.cpp:1612
ProSHADE_settings::changeMapResolution
bool changeMapResolution
Should maps be re-sampled to obtain the required resolution?
Definition: ProSHADE_settings.hpp:51
ProSHADE_internal_tasks::ReportDistancesResults
void ReportDistancesResults(ProSHADE_settings *settings, std::string str1, std::string str2, proshade_double enLevDist, proshade_double trSigmDist, proshade_double rotFunDist)
Simple function for reporting the distances computation results.
Definition: ProSHADE_tasks.cpp:228
ProSHADE_internal_distances::computeTraceSigmaDescriptor
proshade_double computeTraceSigmaDescriptor(ProSHADE_internal_data::ProSHADE_data *obj1, ProSHADE_internal_data::ProSHADE_data *obj2, ProSHADE_settings *settings)
This function computes the trace sigma descriptor value between two objects.
Definition: ProSHADE_distances.cpp:616
ProSHADE_internal_data::ProSHADE_data::getYToPtr
proshade_signed * getYToPtr(void)
This function allows access to the map last position along the Y axis.
Definition: ProSHADE_data.cpp:4009
ProSHADE_internal_data::ProSHADE_data::detectSymmetryInStructure
void detectSymmetryInStructure(ProSHADE_settings *settings, std::vector< proshade_double * > *axes, std::vector< std::vector< proshade_double > > *allCs)
This function runs the symmetry detection algorithms on this structure and saves the results in the s...
Definition: ProSHADE_data.cpp:1801
ProSHADE_internal_data::ProSHADE_data::getXToPtr
proshade_signed * getXToPtr(void)
This function allows access to the map last position along the X axis.
Definition: ProSHADE_data.cpp:3999
ProSHADE_internal_data::ProSHADE_data::getYFromPtr
proshade_signed * getYFromPtr(void)
This function allows access to the map start along the Y axis.
Definition: ProSHADE_data.cpp:3979
ProSHADE_internal_data::ProSHADE_data::mapCOMProcessChangeY
proshade_double mapCOMProcessChangeY
The change in Y axis between the creation of the structure (originalMapYCom) and just before rotation...
Definition: ProSHADE_data.hpp:98
ProSHADE_internal_data::ProSHADE_data::reportSymmetryResults
void reportSymmetryResults(ProSHADE_settings *settings)
This function takes prints the report for symmetry detection.
Definition: ProSHADE_data.cpp:3498
ProSHADE_settings::inputFiles
std::vector< std::string > inputFiles
This vector contains the filenames of all input structure files.
Definition: ProSHADE_settings.hpp:43
ProSHADE_internal_distances::computeEnergyLevelsDescriptor
proshade_double computeEnergyLevelsDescriptor(ProSHADE_internal_data::ProSHADE_data *obj1, ProSHADE_internal_data::ProSHADE_data *obj2, ProSHADE_settings *settings)
This function computes the energy levels descriptor value between two objects.
Definition: ProSHADE_distances.cpp:161
ProSHADE_internal_tasks::checkMapManipulationSettings
void checkMapManipulationSettings(ProSHADE_settings *settings)
The re-boxing settings checks.
Definition: ProSHADE_tasks.cpp:106
ProSHADE_tasks.hpp
This header declares all the taks functions.
ProSHADE_internal_data::ProSHADE_data::originalPdbRotCenZ
proshade_double originalPdbRotCenZ
The centre of rotation as it relates to the original PDB positions (and not the ProSHADE internal map...
Definition: ProSHADE_data.hpp:104
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_internal_tasks::checkSymmetrySettings
void checkSymmetrySettings(ProSHADE_settings *settings)
The symmetry computation settings checks.
Definition: ProSHADE_tasks.cpp:348
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_internal_data::ProSHADE_data::detectSymmetryFromAngleAxisSpace
void detectSymmetryFromAngleAxisSpace(ProSHADE_settings *settings, std::vector< proshade_double * > *axes, std::vector< std::vector< proshade_double > > *allCs)
This function runs the symmetry detection algorithms on this structure using the angle-axis space and...
Definition: ProSHADE_data.cpp:1934
ProSHADE_internal_data::ProSHADE_data::originalPdbRotCenX
proshade_double originalPdbRotCenX
The centre of rotation as it relates to the original PDB positions (and not the ProSHADE internal map...
Definition: ProSHADE_data.hpp:102
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_data::ProSHADE_data::readInStructure
void readInStructure(std::string fName, proshade_unsign inputO, ProSHADE_settings *settings)
This function initialises the basic ProSHADE_data variables and reads in a single structure.
Definition: ProSHADE_data.cpp:501
ProSHADE_internal_misc::deepCopyBoundsSigPtrVector
void deepCopyBoundsSigPtrVector(std::vector< proshade_signed * > *sigPtrVec, proshade_signed *xFrom, proshade_signed *xTo, proshade_signed *yFrom, proshade_signed *yTo, proshade_signed *zFrom, proshade_signed *zTo)
Does a deep copy of a signed int array to a vector of signed int arrays.
Definition: ProSHADE_misc.cpp:347
ProSHADE_internal_tasks::SymmetryDetectionTask
void SymmetryDetectionTask(ProSHADE_settings *settings, std::vector< proshade_double * > *axes, std::vector< std::vector< proshade_double > > *allCs, std::vector< proshade_double > *mapCOMShift)
The symmetry detection task driver function.
Definition: ProSHADE_tasks.cpp:287
ProSHADE_internal_data::ProSHADE_data::getZFromPtr
proshade_signed * getZFromPtr(void)
This function allows access to the map start along the Z axis.
Definition: ProSHADE_data.cpp:3989
ProSHADE_internal_data::ProSHADE_data::getZToPtr
proshade_signed * getZToPtr(void)
This function allows access to the map last position along the Z axis.
Definition: ProSHADE_data.cpp:4019
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