ProSHADE  0.7.6.6 (JUL 2022)
Protein Shape Detection
ProSHADE_internal_tasks Namespace Reference

This namespace contains the main driving functions for each task. More...

Functions

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. More...
 
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. More...
 
void SymmetryDetectionTask (ProSHADE_settings *settings, std::vector< proshade_double > *mapCOMShift, std::string *symT, proshade_unsign *symF, std::vector< proshade_double * > *symA, std::vector< std::vector< proshade_double > > *allCs)
 The symmetry detection task driver function. More...
 
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. More...
 
void SymmetryCentreDetectionTask (ProSHADE_settings *settings, proshade_unsign strIndex=0)
 The task for finding the structure centre based on phase-less symmetry detection. More...
 
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. More...
 
void checkMapManipulationSettings (ProSHADE_settings *settings)
 The re-boxing settings checks. More...
 
void checkDistancesSettings (ProSHADE_settings *settings)
 The distances computation settings checks. More...
 
void checkSymmetrySettings (ProSHADE_settings *settings)
 The symmetry computation settings checks. More...
 
void checkOverlaySettings (ProSHADE_settings *settings)
 The map overlay computation settings checks. More...
 

Detailed Description

This namespace contains the main driving functions for each task.

The ProSHADE_internal_tasks namespace contains the driving functions for all the tasks that can be accomplished by the ProSHADE tool. The user should not need to access this namespace when using the library.

Function Documentation

◆ checkDistancesSettings()

void ProSHADE_internal_tasks::checkDistancesSettings ( ProSHADE_settings settings)

The distances computation settings checks.

This function is called to check the settings object for having all the required information for the distances computation task to proceed.

Parameters
[in]settingsProSHADE_settings object specifying the details of how distances computation should be done.

Definition at line 258 of file ProSHADE_tasks.cpp.

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 }

◆ checkMapManipulationSettings()

void ProSHADE_internal_tasks::checkMapManipulationSettings ( ProSHADE_settings settings)

The re-boxing settings checks.

This function is called to check the settings object for having all the required information for the Re-Boxing task to proceed.

Parameters
[in]settingsProSHADE_settings object specifying the details of how re-boxing should be done.

Definition at line 106 of file ProSHADE_tasks.cpp.

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 }

◆ checkOverlaySettings()

void ProSHADE_internal_tasks::checkOverlaySettings ( ProSHADE_settings settings)

The map overlay computation settings checks.

This function is called to check the settings object for having all the required information for the map overlay task to proceed.

Parameters
[in]settingsProSHADE_settings object specifying the details of how map overlay should be done.

Definition at line 696 of file ProSHADE_tasks.cpp.

697 {
698  //================================================ Are the any structures?
699  if ( settings->inputFiles.size () != 2 )
700  {
701  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." );
702  }
703 
704  //================================================ If centring is on, turn it off and report warning.
705  if ( settings->moveToCOM )
706  {
707  ProSHADE_internal_messages::printWarningMessage ( settings->verbose, "!!! ProSHADE WARNING !!! Map centring was requested, but makes no sense for overlay mode. Turning it off.", "WO00066" );
708  settings->moveToCOM = false;
709  }
710 
711  //================================================ Done
712  return ;
713 
714 }

◆ checkSymmetrySettings()

void ProSHADE_internal_tasks::checkSymmetrySettings ( ProSHADE_settings settings)

The symmetry computation settings checks.

This function is called to check the settings object for having all the required information for the symmetry computation task to proceed.

Parameters
[in]settingsProSHADE_settings object specifying the details of how symmetry detection should be done.

Definition at line 603 of file ProSHADE_tasks.cpp.

604 {
605  //================================================ Are the any structures?
606  if ( settings->inputFiles.size () < 1 )
607  {
608  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." );
609  }
610 
611  //================================================ Is the axis tolerance set properly?
612  if ( settings->axisErrTolerance < 0.0 )
613  {
614  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." );
615  }
616 
617  //================================================ Done
618  return ;
619 
620 }

◆ DistancesComputationTask()

void ProSHADE_internal_tasks::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.

This function is called to proceed with the distances computation task according to the information placed in the settings object passed as the first argument.

Parameters
[in]settingsProSHADE_settings object specifying the details of how distances computation should be done.
[in]enLevsPointer to vector where all energy levels distances are to be saved into.
[in]trSigmPointer to vector where all trace sigma distances are to be saved into.
[in]rotFunPointer to vector where all rotation function distances are to be saved into.

Definition at line 147 of file ProSHADE_tasks.cpp.

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.", settings->messageShift ); }
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.", settings->messageShift ); }
193 
194  proshade_double rotFunDist = 0.0;
195  if ( settings->computeRotationFuncDesc ) { rotFunDist = ProSHADE_internal_distances::computeRotationFunctionDescriptor ( compareAgainst, compareChanging, settings ); }
196  else { ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 1, "Rotation function distance computation not required.", settings->messageShift ); }
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 }

◆ MapManipulationTask()

void ProSHADE_internal_tasks::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.

This function is called to proceed with the map re-boxing task according to the information placed in the settings object passed as the first argument.

Parameters
[in]settingsProSHADE_settings object specifying the details of how re-boxing should be done.
[in]originalBoundsVector to which the original map boundaries of each re-boxed map will be saved into.
[in]reboxedBoundsVector to which the re-boxed map boundaries of each re-boxed map will be saved into.
[in]manipulatedMapsVector to which the map values of each re-boxed map will be saved into.

Definition at line 35 of file ProSHADE_tasks.cpp.

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(), settings->messageShift );
76  if ( settings->reBoxMap ) { reBoxStr->writeMap ( ss.str() ); }
77  else { strToRebox->writeMap ( ss.str() ); }
78  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, "Structure saved.", settings->messageShift );
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 }

◆ MapOverlayTask()

void ProSHADE_internal_tasks::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.

This function is called to run the detect symmetries task according to the information placed in the settings object passed as the first argument.

Parameters
[in]settingsProSHADE_settings object specifying the details of how distances computation should be done.
[in]rotationCentrePointer to vector for saving the position of the centre of rotation about which the rotation is to be done.
[in]eulerAnglesPointer to vector where the three Euler angles will be saved into.
[in]finalTranslationPointer to a vector where the translation required to move structure from origin to optimal overlay with static structure will be saved into.

Definition at line 632 of file ProSHADE_tasks.cpp.

633 {
634  //================================================ Check the settings are complete and meaningful
635  checkOverlaySettings ( settings );
636 
637  //================================================ Initialise variables
638  proshade_double eulA, eulB, eulG, trsX, trsY, trsZ;
639 
640  //================================================ Create the data objects initially (this time without phase)
643 
644  //================================================ First, run without phase and find best rotation angles
645  settings->usePhase = false;
646  ProSHADE_internal_overlay::getOptimalRotation ( settings, staticStructure, movingStructure, &eulA, &eulB, &eulG );
647 
648  //================================================ Release memory
649  delete staticStructure;
650  delete movingStructure;
651 
652  //================================================ Create the data objects again (this time with phase)
653  staticStructure = new ProSHADE_internal_data::ProSHADE_data ( );
654  movingStructure = new ProSHADE_internal_data::ProSHADE_data ( );
655 
656  //================================================ Now, run with phase and find optimal translation
657  settings->usePhase = true;
658  settings->changeMapResolution = true;
659  ProSHADE_internal_overlay::getOptimalTranslation ( settings, staticStructure, movingStructure, &trsX, &trsY, &trsZ, eulA, eulB, eulG );
660 
661  //================================================ Compute the proper translations using the translation function output
662  ProSHADE_internal_misc::addToDoubleVector ( rotationCentre, movingStructure->originalPdbRotCenX );
663  ProSHADE_internal_misc::addToDoubleVector ( rotationCentre, movingStructure->originalPdbRotCenY );
664  ProSHADE_internal_misc::addToDoubleVector ( rotationCentre, movingStructure->originalPdbRotCenZ );
665  ProSHADE_internal_misc::addToDoubleVector ( finalTranslation, movingStructure->originalPdbTransX );
666  ProSHADE_internal_misc::addToDoubleVector ( finalTranslation, movingStructure->originalPdbTransY );
667  ProSHADE_internal_misc::addToDoubleVector ( finalTranslation, movingStructure->originalPdbTransZ );
668 
669  //================================================ Write out everything
670  movingStructure->writeOutOverlayFiles ( settings, eulA, eulB, eulG, rotationCentre, finalTranslation );
671 
672  //================================================ Save the rotation and rest of translations
673  ProSHADE_internal_misc::addToDoubleVector ( eulerAngles, eulA );
674  ProSHADE_internal_misc::addToDoubleVector ( eulerAngles, eulB );
675  ProSHADE_internal_misc::addToDoubleVector ( eulerAngles, eulG );
676 
677  //================================================ Report results to user
678  movingStructure->reportOverlayResults ( settings, rotationCentre, eulerAngles, finalTranslation );
679 
680  //================================================ Release memory
681  delete staticStructure;
682  delete movingStructure;
683 
684  //================================================ Done
685  return ;
686 
687 }

◆ ReportDistancesResults()

void ProSHADE_internal_tasks::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.

Parameters
[in]settingsProSHADE_settings object specifying the details of how distances computation should be done.
[in]str1The name of the structure to which all other structures are to be compared to.
[in]str2The name of the structure which is compared to str1.
[in]enLevDistThe value of the energy levels descriptor for the two structures.
[in]trSimDistThe value of the trace sigma descriptor for the two structures.
[in]rotFunDistThe value of the roation function descriptor for the two structures.

Definition at line 228 of file ProSHADE_tasks.cpp.

229 {
230  std::stringstream hlpSS;
231  hlpSS << "Distances between " << str1 << " and " << str2;
232  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 0, hlpSS.str(), settings->messageShift );
233 
234  std::stringstream hlpSSE;
235  hlpSSE << "Energy levels distance : " << enLevDist;
236  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 0, hlpSSE.str(), settings->messageShift );
237 
238  std::stringstream hlpSSS;
239  hlpSSS << "Trace sigma distance : " << trSigmDist;
240  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 0, hlpSSS.str(), settings->messageShift );
241 
242  std::stringstream hlpSSR;
243  hlpSSR << "Rotation function distance: " << rotFunDist;
244  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 0, hlpSSR.str(), settings->messageShift );
245 
246  //================================================ Done
247  return ;
248 
249 }

◆ SymmetryCentreDetectionTask()

void ProSHADE_internal_tasks::SymmetryCentreDetectionTask ( ProSHADE_settings settings,
proshade_unsign  strIndex = 0 
)

The task for finding the structure centre based on phase-less symmetry detection.

This function is called to compute the symmetry of the phase-less map so that (in case there is any) it could then find the centre of rotation and thus the centre of the structure.

Parameters
[in]settingsProSHADE_settings object specifying the details of how symmetry centre detection should be done.
[in]strIndexThe index of the structure to be read from the structure list available in the settings object.

Definition at line 374 of file ProSHADE_tasks.cpp.

375 {
376  //================================================ Keep original settings for the phased reading
377  ProSHADE_settings* tmpSettings = new ProSHADE_settings ( settings );
378 
379  //================================================ Enforce the necessary settings
380  tmpSettings->usePhase = false;
381  tmpSettings->requestedSymmetryType = "onlyCandD";
382  tmpSettings->moveToCOM = false;
383  tmpSettings->addExtraSpace = tmpSettings->addExtraSpace * 5.0f;
384  settings->moveToCOM = false;
385 
386  //================================================ Read in the structure and find all symmetries without using phase information
388  symStr->readInStructure ( tmpSettings->inputFiles.at(strIndex), strIndex, tmpSettings );
389  symStr->processInternalMap ( tmpSettings );
390  symStr->mapToSpheres ( tmpSettings );
391  symStr->computeSphericalHarmonics ( tmpSettings );
392  symStr->computeRotationFunction ( tmpSettings );
393  symStr->detectSymmetryFromAngleAxisSpace ( tmpSettings );
394 
395  //================================================ Find reliable symmetries in the Patterson map
396  std::vector< proshade_unsign > relSym = ProSHADE_internal_symmetry::findReliableUnphasedSymmetries ( symStr->getCyclicAxes(), symStr->getDihedralAxes(), tmpSettings->verbose, settings->messageShift, tmpSettings->axisErrTolerance );
397  std::vector< std::vector< proshade_double > > allCs;
398  std::vector< proshade_double > hlpVec;
399  for ( size_t it1 = 0; it1 < symStr->getCyclicAxes()->size(); it1++ )
400  {
401  hlpVec.clear ( );
402  for ( size_t it2 = 0; it2 < 7; it2++ ) { ProSHADE_internal_misc::addToDoubleVector ( &hlpVec, symStr->getCyclicAxes()->at(it1)[it2] ); }
404  }
405 
406  //================================================ If no symmetries are found, inform the user
407  if ( relSym.size() == 0 )
408  {
409  ProSHADE_internal_messages::printWarningMessage ( tmpSettings->verbose, "!!! ProSHADE WARNING !!! Failed to find symmetry in Patterson map. Map rotation centre detection cannot be done without a symmetry, returning vector with [Inf, Inf, Inf].", "WS00071" );
410  settings->centrePosition.at(0) = std::numeric_limits< proshade_double >::infinity();
411  settings->centrePosition.at(1) = std::numeric_limits< proshade_double >::infinity();
412  settings->centrePosition.at(2) = std::numeric_limits< proshade_double >::infinity();
413  return ;
414  }
415  else
416  {
417  std::stringstream ssHlp;
418  std::stringstream ssHlp2;
419  if ( relSym.size() == 1 )
420  {
421  ssHlp << "Decided that there is one reliable axis forming C" << allCs.at(relSym.at(0))[0];
422  ssHlp2 << "The axis is: [" << allCs.at(relSym.at(0))[1] << " , " << allCs.at(relSym.at(0))[2] << " , " << allCs.at(relSym.at(0))[3] << "] with angle " << allCs.at(relSym.at(0))[4] << " peak height " << allCs.at(relSym.at(0))[5] << " and average FSC of " << allCs.at(relSym.at(0))[6];
423  }
424  else
425  {
426  ssHlp << "Decided that there are two reliable axes forming D" << std::max ( allCs.at(relSym.at(0))[0], allCs.at(relSym.at(1))[0] );
427  ssHlp2 << "The first axis is: [" << allCs.at(relSym.at(0))[1] << " , " << allCs.at(relSym.at(0))[2] << " , " << allCs.at(relSym.at(0))[3] << "] with angle " << allCs.at(relSym.at(0))[4] << " peak height " << allCs.at(relSym.at(0))[5] << " and average FSC of " << allCs.at(relSym.at(0))[6] << std::endl;
428  for ( size_t it = 0; it < static_cast< size_t > ( settings->verbose + settings->messageShift ); it++ ) { ssHlp2 << " ..."; }
429  ssHlp2 << " and the second axis is: [" << allCs.at(relSym.at(1))[1] << " , " << allCs.at(relSym.at(1))[2] << " , " << allCs.at(relSym.at(1))[3] << "] with angle " << allCs.at(relSym.at(1))[4] << " peak height " << allCs.at(relSym.at(1))[5] << " and average FSC of " << allCs.at(relSym.at(1))[6];
430  }
431  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, ssHlp.str(), settings->messageShift );
432  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 5, ssHlp2.str(), settings->messageShift );
433 
434 
435  }
436 
437  //================================================ Optimise the orthogonal pair, if there is one
438  if ( relSym.size() == 2 )
439  {
440  //============================================ Optimise the orthogonal pair
441  ProSHADE_internal_symmetry::optimiseDGroupAngleFromAxesHeights ( &allCs, relSym, symStr, tmpSettings );
442  }
443 
444  //================================================ Re-read the map, this time with phases
445  delete symStr;
446  symStr = new ProSHADE_internal_data::ProSHADE_data ( );
447  symStr->readInStructure ( settings->inputFiles.at(strIndex), strIndex, settings );
448  symStr->processInternalMap ( settings );
449 
450  //================================================ Allocate the Fourier transforms related memory
451  fftw_complex *origMap = nullptr, *origCoeffs = nullptr, *rotMapComplex = nullptr, *rotCoeffs = nullptr, *trFunc = nullptr, *trFuncCoeffs = nullptr;
452  fftw_plan planForwardFourier, planForwardFourierRot, planReverseFourierComb;
453  ProSHADE_internal_symmetry::allocateCentreOfMapFourierTransforms ( symStr->getXDim(), symStr->getYDim(), symStr->getZDim(), origMap, origCoeffs, rotMapComplex, rotCoeffs, trFunc, trFuncCoeffs, &planForwardFourier, &planForwardFourierRot, &planReverseFourierComb );
454 
455  //================================================ Compute Fourier for the original map
456  for ( proshade_unsign it = 0; it < static_cast< proshade_unsign > ( symStr->getXDim() * symStr->getYDim() * symStr->getZDim() ); it++ ) { origMap[it][0] = symStr->getMapValue( it ); origMap[it][1] = 0.0; }
457  fftw_execute ( planForwardFourier );
458 
459  //================================================ If single C was found
460  if ( relSym.size() == 1 )
461  {
462  //============================================ Initialise local variables
463  proshade_double xMapCOM = 0.0, yMapCOM = 0.0, zMapCOM = 0.0;
464  std::vector< proshade_unsign > axLst;
465  std::vector< std::vector< proshade_double > > symElems;
466 
467  //============================================ Find the point
468  ProSHADE_internal_misc::addToUnsignVector ( &axLst, static_cast< proshade_unsign > ( relSym.at(0) ) );
469  symElems = symStr->getAllGroupElements ( &allCs, axLst, "C", tmpSettings->axisErrTolerance );
470  std::vector< proshade_double > pointPos = ProSHADE_internal_symmetry::findPointFromTranslations ( symStr,
471  symElems,
472  origCoeffs, rotMapComplex,
473  rotCoeffs, planForwardFourierRot,
474  trFuncCoeffs, trFunc,
475  planReverseFourierComb );
476 
477  //============================================ Find COM in Angstroms in visualisation space
478  ProSHADE_internal_mapManip::findMAPCOMValues ( symStr->internalMap, &xMapCOM, &yMapCOM, &zMapCOM, symStr->xDimSize, symStr->yDimSize, symStr->zDimSize, symStr->xFrom, symStr->xTo, symStr->yFrom, symStr->yTo, symStr->zFrom, symStr->zTo, settings->removeNegativeDensity );
479 
480  //============================================ Determine box centre in indices
481  proshade_double xBoxCentre = ( ( static_cast< proshade_double > ( symStr->xTo ) - static_cast< proshade_double > ( symStr->xFrom ) ) / 2.0 ) + static_cast< proshade_double > ( symStr->xFrom );
482  proshade_double yBoxCentre = ( ( static_cast< proshade_double > ( symStr->yTo ) - static_cast< proshade_double > ( symStr->yFrom ) ) / 2.0 ) + static_cast< proshade_double > ( symStr->yFrom );
483  proshade_double zBoxCentre = ( ( static_cast< proshade_double > ( symStr->zTo ) - static_cast< proshade_double > ( symStr->zFrom ) ) / 2.0 ) + static_cast< proshade_double > ( symStr->zFrom );
484 
485  //============================================ Determine distance from COM in indices to box centre in indices
486  proshade_double xCOMFromBoxCen = xBoxCentre - ( xMapCOM / ( static_cast< proshade_double > ( symStr->xDimSize ) / static_cast< proshade_double > ( symStr->xDimIndices ) ) );
487  proshade_double yCOMFromBoxCen = yBoxCentre - ( yMapCOM / ( static_cast< proshade_double > ( symStr->yDimSize ) / static_cast< proshade_double > ( symStr->yDimIndices ) ) );
488  proshade_double zCOMFromBoxCen = zBoxCentre - ( zMapCOM / ( static_cast< proshade_double > ( symStr->zDimSize ) / static_cast< proshade_double > ( symStr->zDimIndices ) ) );
489 
490  //============================================ Determine the coefficient of mapping of the COM point to the line
491  proshade_double alpha1 = ProSHADE_internal_maths::computeDotProduct ( pointPos.at(0) - xCOMFromBoxCen,
492  pointPos.at(1) - yCOMFromBoxCen,
493  pointPos.at(2) - zCOMFromBoxCen,
494  allCs.at(relSym.at(0))[1],
495  allCs.at(relSym.at(0))[2],
496  allCs.at(relSym.at(0))[3] ) /
497  ProSHADE_internal_maths::computeDotProduct ( allCs.at(relSym.at(0))[1],
498  allCs.at(relSym.at(0))[2],
499  allCs.at(relSym.at(0))[3],
500  allCs.at(relSym.at(0))[1],
501  allCs.at(relSym.at(0))[2],
502  allCs.at(relSym.at(0))[3] );
503 
504  //============================================ Determine point on the axis closest to COM
505  settings->centrePosition.at(0) = pointPos.at(0) + ( alpha1 * allCs.at(relSym.at(0))[1] );
506  settings->centrePosition.at(1) = pointPos.at(1) + ( alpha1 * allCs.at(relSym.at(0))[2] );
507  settings->centrePosition.at(2) = pointPos.at(2) + ( alpha1 * allCs.at(relSym.at(0))[3] );
508  }
509  //================================================ If D was found
510  else
511  {
512  //============================================ Initialise local variables
513  std::vector< proshade_unsign > axLst;
514  std::vector< std::vector< proshade_double > > symElems;
515 
516  //============================================ Find the first point
517  ProSHADE_internal_misc::addToUnsignVector ( &axLst, static_cast< proshade_unsign > ( relSym.at(0) ) );
518  symElems = symStr->getAllGroupElements ( &allCs, axLst, "C", tmpSettings->axisErrTolerance );
519  std::vector< proshade_double > point1Pos = ProSHADE_internal_symmetry::findPointFromTranslations ( symStr,
520  symElems,
521  origCoeffs, rotMapComplex,
522  rotCoeffs, planForwardFourierRot,
523  trFuncCoeffs, trFunc,
524  planReverseFourierComb );
525 
526  //============================================ Find the second point
527  axLst.at(0) = static_cast< proshade_unsign > ( relSym.at(1) );
528  symElems = symStr->getAllGroupElements ( &allCs, axLst, "C", tmpSettings->axisErrTolerance );
529  std::vector< proshade_double > point2Pos = ProSHADE_internal_symmetry::findPointFromTranslations ( symStr,
530  symElems,
531  origCoeffs, rotMapComplex,
532  rotCoeffs, planForwardFourierRot,
533  trFuncCoeffs, trFunc,
534  planReverseFourierComb );
535 
536  //============================================ Compute the tangents
537  proshade_double* tangentToAxes = ProSHADE_internal_maths::computeCrossProduct ( allCs.at(relSym.at(0))[1], allCs.at(relSym.at(0))[2], allCs.at(relSym.at(0))[3],
538  allCs.at(relSym.at(1))[1], allCs.at(relSym.at(1))[2], allCs.at(relSym.at(1))[3] );
539 
540  proshade_double* correctedSecondAxis = ProSHADE_internal_maths::computeCrossProduct ( allCs.at(relSym.at(0))[1], allCs.at(relSym.at(0))[2], allCs.at(relSym.at(0))[3],
541  tangentToAxes[0], tangentToAxes[1], tangentToAxes[2] );
542 
543  proshade_double* correctedFirstAxis = ProSHADE_internal_maths::computeCrossProduct ( allCs.at(relSym.at(1))[1], allCs.at(relSym.at(1))[2], allCs.at(relSym.at(1))[3],
544  tangentToAxes[0], tangentToAxes[1], tangentToAxes[2] );
545 
546  //============================================ Compute mappings
547  proshade_double alpha1 = ProSHADE_internal_maths::computeDotProduct ( point2Pos.at(0) - point1Pos.at(0),
548  point2Pos.at(1) - point1Pos.at(1),
549  point2Pos.at(2) - point1Pos.at(2),
550  correctedFirstAxis[0], correctedFirstAxis[1], correctedFirstAxis[2] ) /
551  ProSHADE_internal_maths::computeDotProduct ( allCs.at(relSym.at(0))[1], allCs.at(relSym.at(0))[2], allCs.at(relSym.at(0))[3],
552  correctedFirstAxis[0], correctedFirstAxis[1], correctedFirstAxis[2] );
553  proshade_double alpha2 = ProSHADE_internal_maths::computeDotProduct ( point1Pos.at(0) - point2Pos.at(0),
554  point1Pos.at(1) - point2Pos.at(1),
555  point1Pos.at(2) - point2Pos.at(2),
556  correctedSecondAxis[0], correctedSecondAxis[1], correctedSecondAxis[2] ) /
557  ProSHADE_internal_maths::computeDotProduct ( allCs.at(relSym.at(1))[1], allCs.at(relSym.at(1))[2], allCs.at(relSym.at(1))[3],
558  correctedSecondAxis[0], correctedSecondAxis[1], correctedSecondAxis[2] );
559 
560  //============================================ Find the intersect point (averaged)
561  settings->centrePosition.at(0) = ( ( point1Pos.at(0) + ( alpha1 * allCs.at(relSym.at(0))[1] ) ) + ( point2Pos.at(0) + ( alpha2 * allCs.at(relSym.at(1))[1] ) ) ) / 2.0;
562  settings->centrePosition.at(1) = ( ( point1Pos.at(1) + ( alpha1 * allCs.at(relSym.at(0))[2] ) ) + ( point2Pos.at(1) + ( alpha2 * allCs.at(relSym.at(1))[2] ) ) ) / 2.0;
563  settings->centrePosition.at(2) = ( ( point1Pos.at(2) + ( alpha1 * allCs.at(relSym.at(0))[3] ) ) + ( point2Pos.at(2) + ( alpha2 * allCs.at(relSym.at(1))[3] ) ) ) / 2.0;
564 
565  //============================================ Release memory
566  delete[] tangentToAxes;
567  delete[] correctedSecondAxis;
568  delete[] correctedFirstAxis;
569 
570  }
571 
572  //================================================ Release the Fourier transforms related memory
573  ProSHADE_internal_symmetry::releaseCentreOfMapFourierTransforms ( origMap, origCoeffs, rotMapComplex, rotCoeffs, trFunc, trFuncCoeffs, planForwardFourier, planForwardFourierRot, planReverseFourierComb );
574 
575  //== Release optimisation memory
576 // delete[] trsOptMap;
577 // delete[] trsOptCoeffs;
578 // fftw_destroy_plan ( planForwardOptimisation );
579 
580  //================================================ Release memory after FSC computation
581 // delete[] FSCmapData;
582 // delete[] FSCorigCoeffs;
583 // delete[] FSCfCoeffs;
584 // fftw_destroy_plan ( FSCplanForwardFourier );
585 // delete[] binIndexing;
586 // for (size_t binIt = 0; binIt < static_cast< size_t > ( noBins ); binIt++ ) { delete[] binDataFSC[binIt]; }
587 // delete[] binDataFSC;
588 // delete[] binCounts;
589 // delete[] fscByBin;
590 
591  //================================================ Done
592  return ;
593 
594 }

◆ SymmetryDetectionTask()

void ProSHADE_internal_tasks::SymmetryDetectionTask ( ProSHADE_settings settings,
std::vector< proshade_double > *  mapCOMShift,
std::string *  symT,
proshade_unsign *  symF,
std::vector< proshade_double * > *  symA,
std::vector< std::vector< proshade_double > > *  allCs 
)

The symmetry detection task driver function.

This function is called to run the detect symmetries task according to the information placed in the settings object passed as the first argument.

Parameters
[in]settingsProSHADE_settings object specifying the details of how distances computation should be done.
[in]mapCOMShiftA pointer to a vector containing the distance from the centre of the map to the point about which the symmetry detection was done.
[in]symTPointer to string that will hold the determined symmetry type.
[in]symFPointer to unsigned int that will hold the determined symmetry fold.
[in]symAPointer to a vector of doubles that will hold the determined symmetry axes in the proshade format.
[in]allCsA pointer to a vector of vectors containing all the detected C axes.

Definition at line 290 of file ProSHADE_tasks.cpp.

291 {
292  //================================================ Check the settings are complete and meaningful
293  checkSymmetrySettings ( settings );
294 
295  //================================================ Now, for each other structure
296  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( settings->inputFiles.size() ); iter++ )
297  {
298  //============================================ Create a data object
300 
301  //============================================ Read in the compared structure
302  symmetryStructure->readInStructure ( settings->inputFiles.at(iter), iter, settings );
303 
304  //============================================ Assume symmetry centre at the box centre, or find it out using Patterson map?
305  if ( settings->findSymCentre )
306  {
307  //======================================== Report progress
308  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 1, "Starting symmetry centre detection procedure.", settings->messageShift );
309 
310  //======================================== Make a local copy of settings (to avoid centre detection settings things for the symmetry detection which will follow)
311  ProSHADE_settings* rotCenSettings = new ProSHADE_settings ( settings );
312  rotCenSettings->messageShift = 1;
313 
314  //======================================== Run the detection
315  SymmetryCentreDetectionTask ( rotCenSettings, iter );
316 
317  //======================================== Save the results
318  settings->centrePosition.at(0) = rotCenSettings->centrePosition.at(0);
319  settings->centrePosition.at(1) = rotCenSettings->centrePosition.at(1);
320  settings->centrePosition.at(2) = rotCenSettings->centrePosition.at(2);
321 
322  //======================================== Report progress
323  std::stringstream ss;
324  ss << "Detected symmetry centre at " << settings->centrePosition.at(0) << " ; " << settings->centrePosition.at(1) << " ; " << settings->centrePosition.at(2);
325  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, ss.str(), settings->messageShift );
326  }
327 
328  //============================================ Internal data processing (COM, norm, mask, extra space)
329  symmetryStructure->processInternalMap ( settings );
330 
331  //============================================ Map to sphere
332  symmetryStructure->mapToSpheres ( settings );
333 
334  //============================================ Get spherical harmonics
335  symmetryStructure->computeSphericalHarmonics ( settings );
336 
337  //============================================ Compute auto-rotation map
338  symmetryStructure->computeRotationFunction ( settings );
339 
340  //======================================== Detect point groups in the angle-axis space
341  symmetryStructure->detectSymmetryFromAngleAxisSpace ( settings );
342 
343  //============================================ Report results
344  symmetryStructure->reportSymmetryResultsList ( settings );
345 
346  //============================================ Save symmetry results
347  *symT = symmetryStructure->getRecommendedSymmetryType ( );
348  *symF = symmetryStructure->getRecommendedSymmetryFold ( );
349  for ( size_t aIt = 0; aIt < symmetryStructure->recommendedSymmetryValues.size(); aIt++ ) { ProSHADE_internal_misc::deepCopyAxisToDblPtrVector ( symA, &symmetryStructure->recommendedSymmetryValues.at( aIt )[0] ); }
350  for ( size_t aIt = 0; aIt < symmetryStructure->cyclicSymmetries.size(); aIt++ ) { std::vector< proshade_double > hlpVec; for ( size_t vIt = 0; vIt < 7; vIt++ ) { ProSHADE_internal_misc::addToDoubleVector ( &hlpVec, symmetryStructure->cyclicSymmetries.at(aIt)[vIt] ); } ProSHADE_internal_misc::addToDoubleVectorVector ( allCs, hlpVec ); }
351 
352  //============================================ Save internal map shift to run object
353  ProSHADE_internal_misc::addToDoubleVector ( mapCOMShift, symmetryStructure->mapCOMProcessChangeX );
354  ProSHADE_internal_misc::addToDoubleVector ( mapCOMShift, symmetryStructure->mapCOMProcessChangeY );
355  ProSHADE_internal_misc::addToDoubleVector ( mapCOMShift, symmetryStructure->mapCOMProcessChangeZ );
356 
357  //============================================ Release memory
358  delete symmetryStructure;
359  }
360 
361  //================================================ Done
362  return ;
363 
364 }
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:1895
ProSHADE_internal_data::ProSHADE_data::cyclicSymmetries
std::vector< proshade_double * > cyclicSymmetries
This is where the detected cyclic ("C") symmetries will be kept.
Definition: ProSHADE_data.hpp:144
ProSHADE_internal_data::ProSHADE_data::zFrom
proshade_signed zFrom
This is the starting index along the z axis.
Definition: ProSHADE_data.hpp:112
ProSHADE_internal_data::ProSHADE_data::xDimIndices
proshade_unsign xDimIndices
This is the size of the map cell x dimension in indices.
Definition: ProSHADE_data.hpp:65
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:119
ProSHADE_settings::findSymCentre
bool findSymCentre
Should phase-less map be used to determine centre of symmetry?
Definition: ProSHADE_settings.hpp:130
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:120
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:1410
ProSHADE_internal_io::isFilePDB
bool isFilePDB(std::string fName)
Function determining if the input data type is PDB.
Definition: ProSHADE_io.cpp:38
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:1352
ProSHADE_internal_maths::computeCrossProduct
proshade_double * computeCrossProduct(proshade_double *x1, proshade_double *y1, proshade_double *z1, proshade_double *x2, proshade_double *y2, proshade_double *z2)
Simple 3D vector cross product computation.
Definition: ProSHADE_maths.cpp:1817
ProSHADE_internal_tasks::checkDistancesSettings
void checkDistancesSettings(ProSHADE_settings *settings)
The distances computation settings checks.
Definition: ProSHADE_tasks.cpp:258
ProSHADE_internal_data::ProSHADE_data::getDihedralAxes
std::vector< std::vector< proshade_double * > > * getDihedralAxes(void)
This function allows access to the list of detected dihedral axes.
Definition: ProSHADE_data.cpp:4271
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:114
ProSHADE_internal_symmetry::findPointFromTranslations
std::vector< proshade_double > findPointFromTranslations(ProSHADE_internal_data::ProSHADE_data *symStr, std::vector< std::vector< proshade_double > > symElems, fftw_complex *origCoeffs, fftw_complex *rotMapComplex, fftw_complex *rotCoeffs, fftw_plan planForwardFourierRot, fftw_complex *trFuncCoeffs, fftw_complex *trFunc, fftw_plan planReverseFourierComb)
This function computes the average of optimal translations for a cyclic point group.
Definition: ProSHADE_symmetry.cpp:3985
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:135
ProSHADE_internal_symmetry::allocateCentreOfMapFourierTransforms
void allocateCentreOfMapFourierTransforms(proshade_unsign xDim, proshade_unsign yDim, proshade_unsign zDim, fftw_complex *&origMap, fftw_complex *&origCoeffs, fftw_complex *&rotMapComplex, fftw_complex *&rotCoeffs, fftw_complex *&trFunc, fftw_complex *&trFuncCoeffs, fftw_plan *planForwardFourier, fftw_plan *planForwardFourierRot, fftw_plan *planReverseFourierComb)
This function allocates the required memory for the Fourier transforms required to find the centre of...
Definition: ProSHADE_symmetry.cpp:3841
ProSHADE_internal_distances::computeRotationFunctionDescriptor
proshade_double computeRotationFunctionDescriptor(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:939
ProSHADE_internal_data::ProSHADE_data
This class contains all inputed and derived data for a single structure.
Definition: ProSHADE_data.hpp:49
ProSHADE_settings::removeNegativeDensity
bool removeNegativeDensity
Should the negative density be removed from input files?
Definition: ProSHADE_settings.hpp:47
ProSHADE_internal_data::ProSHADE_data::zDimSize
proshade_single zDimSize
This is the size of the map cell z dimension in Angstroms.
Definition: ProSHADE_data.hpp:61
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:696
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:64
ProSHADE_internal_data::ProSHADE_data::getRecommendedSymmetryFold
proshade_unsign getRecommendedSymmetryFold(void)
This function simply returns the detected recommended symmetry fold.
Definition: ProSHADE_data.cpp:4609
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:74
ProSHADE_internal_symmetry::optimiseDGroupAngleFromAxesHeights
void optimiseDGroupAngleFromAxesHeights(std::vector< std::vector< proshade_double > > *ret, ProSHADE_internal_data::ProSHADE_data *dataObj, ProSHADE_settings *settings)
This function takes two axes with almost dihedral angle and optimises their relative positions as wel...
Definition: ProSHADE_symmetry.cpp:3269
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:102
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:4104
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:152
ProSHADE_internal_data::ProSHADE_data::getAllGroupElements
std::vector< std::vector< proshade_double > > getAllGroupElements(std::vector< proshade_unsign > axesList, std::string groupType="", proshade_double matrixTolerance=0.05)
This function returns the group elements as rotation matrices of any defined point group.
Definition: ProSHADE_data.cpp:3122
ProSHADE_internal_data::ProSHADE_data::yDimIndices
proshade_unsign yDimIndices
This is the size of the map cell y dimension in indices.
Definition: ProSHADE_data.hpp:66
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:3431
ProSHADE_settings::messageShift
proshade_signed messageShift
This value allows shifting the messages to create more readable log for sub-processes.
Definition: ProSHADE_settings.hpp:153
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::getRecommendedSymmetryType
std::string getRecommendedSymmetryType(void)
This function simply returns the detected recommended symmetry type.
Definition: ProSHADE_data.cpp:4600
ProSHADE_internal_data::ProSHADE_data::getXDim
proshade_unsign getXDim(void)
This function allows access to the map size in indices along the X axis.
Definition: ProSHADE_data.cpp:4074
ProSHADE_internal_data::ProSHADE_data::getCyclicAxes
std::vector< proshade_double * > * getCyclicAxes(void)
This function allows access to the list of detected cyclic axes.
Definition: ProSHADE_data.cpp:4239
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::zTo
proshade_signed zTo
This is the final index along the z axis.
Definition: ProSHADE_data.hpp:115
ProSHADE_internal_data::ProSHADE_data::readInStructure
void readInStructure(std::string fName, proshade_unsign inputO, ProSHADE_settings *settings, proshade_double *maskArr=nullptr, proshade_unsign maskXDim=0, proshade_unsign maskYDim=0, proshade_unsign maskZDim=0, proshade_double *weightsArr=nullptr, proshade_unsign weigXDim=0, proshade_unsign weigYDim=0, proshade_unsign weigZDim=0)
This function initialises the basic ProSHADE_data variables and reads in a single structure.
Definition: ProSHADE_data.cpp:525
ProSHADE_internal_data::ProSHADE_data::xTo
proshade_signed xTo
This is the final index along the x axis.
Definition: ProSHADE_data.hpp:113
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:1851
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_data::ProSHADE_data::yTo
proshade_signed yTo
This is the final index along the y axis.
Definition: ProSHADE_data.hpp:114
ProSHADE_internal_data::ProSHADE_data::xDimSize
proshade_single xDimSize
This is the size of the map cell x dimension in Angstroms.
Definition: ProSHADE_data.hpp:59
ProSHADE_internal_maths::computeDotProduct
proshade_double computeDotProduct(proshade_double *x1, proshade_double *y1, proshade_double *z1, proshade_double *x2, proshade_double *y2, proshade_double *z2)
Simple 3D vector dot product computation.
Definition: ProSHADE_maths.cpp:1785
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:433
ProSHADE_internal_data::ProSHADE_data::internalMap
proshade_double * internalMap
The internal map data representation, which may be amended as the run progresses.
Definition: ProSHADE_data.hpp:56
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:4743
ProSHADE_internal_data::ProSHADE_data::xFrom
proshade_signed xFrom
This is the starting index along the x axis.
Definition: ProSHADE_data.hpp:110
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:916
ProSHADE_settings::reBoxMap
bool reBoxMap
This switch decides whether re-boxing is needed.
Definition: ProSHADE_settings.hpp:96
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:123
ProSHADE_settings::addExtraSpace
proshade_single addExtraSpace
If this value is non-zero, this many angstroms of empty space will be added to the internal map.
Definition: ProSHADE_settings.hpp:107
ProSHADE_internal_data::ProSHADE_data::getYDim
proshade_unsign getYDim(void)
This function allows access to the map size in indices along the Y axis.
Definition: ProSHADE_data.cpp:4084
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:4709
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:1739
ProSHADE_internal_mapManip::findMAPCOMValues
void findMAPCOMValues(proshade_double *map, proshade_double *xCom, proshade_double *yCom, proshade_double *zCom, proshade_single xAngs, proshade_single yAngs, proshade_single zAngs, proshade_signed xFrom, proshade_signed xTo, proshade_signed yFrom, proshade_signed yTo, proshade_signed zFrom, proshade_signed zTo, bool removeNegDens)
This function finds the Centre of Mass for a map.
Definition: ProSHADE_mapManip.cpp:243
ProSHADE_internal_symmetry::findReliableUnphasedSymmetries
std::vector< proshade_unsign > findReliableUnphasedSymmetries(std::vector< proshade_double * > *allCs, std::vector< std::vector< proshade_double * > > *allDs, proshade_signed verbose, proshade_signed messageShift, proshade_double tolerance)
This function checks the list of detected axes (presumably from phaseless symmetry detection) and ret...
Definition: ProSHADE_symmetry.cpp:3735
ProSHADE_internal_data::ProSHADE_data::getZDim
proshade_unsign getZDim(void)
This function allows access to the map size in indices along the Z axis.
Definition: ProSHADE_data.cpp:4094
ProSHADE_settings::centrePosition
std::vector< proshade_double > centrePosition
The position of the centre of the map in "real space" co-ordinates.
Definition: ProSHADE_settings.hpp:145
ProSHADE_settings::changeMapResolution
bool changeMapResolution
Should maps be re-sampled to obtain the required resolution?
Definition: ProSHADE_settings.hpp:51
ProSHADE_internal_data::ProSHADE_data::recommendedSymmetryValues
std::vector< proshade_double * > recommendedSymmetryValues
The axes and other info of the recommended symmetry for the structure.
Definition: ProSHADE_data.hpp:149
ProSHADE_internal_data::ProSHADE_data::getMapValue
proshade_double getMapValue(proshade_unsign pos)
This function returns the internal map representation value of a particular array position.
Definition: ProSHADE_data.cpp:3738
ProSHADE_internal_symmetry::releaseCentreOfMapFourierTransforms
void releaseCentreOfMapFourierTransforms(fftw_complex *origMap, fftw_complex *origCoeffs, fftw_complex *rotMapComplex, fftw_complex *rotCoeffs, fftw_complex *trFunc, fftw_complex *trFuncCoeffs, fftw_plan planForwardFourier, fftw_plan planForwardFourierRot, fftw_plan planReverseFourierComb)
This function releases the allocated memory for the Fourier transforms used to find the centre of the...
Definition: ProSHADE_symmetry.cpp:3881
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:648
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:4144
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:103
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:4134
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:4114
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_misc::addToDoubleVectorVector
void addToDoubleVectorVector(std::vector< std::vector< proshade_double > > *vecToAddTo, std::vector< proshade_double > elementToAdd)
Adds the element to the vector of vectors.
Definition: ProSHADE_misc.cpp:233
ProSHADE_internal_data::ProSHADE_data::yFrom
proshade_signed yFrom
This is the starting index along the y axis.
Definition: ProSHADE_data.hpp:111
ProSHADE_internal_misc::addToUnsignVector
void addToUnsignVector(std::vector< proshade_unsign > *vecToAddTo, proshade_unsign elementToAdd)
Adds the element to the vector.
Definition: ProSHADE_misc.cpp:99
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_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:132
ProSHADE_internal_tasks::checkSymmetrySettings
void checkSymmetrySettings(ProSHADE_settings *settings)
The symmetry computation settings checks.
Definition: ProSHADE_tasks.cpp:603
ProSHADE_internal_tasks::SymmetryCentreDetectionTask
void SymmetryCentreDetectionTask(ProSHADE_settings *settings, proshade_unsign strIndex=0)
The task for finding the structure centre based on phase-less symmetry detection.
Definition: ProSHADE_tasks.cpp:374
ProSHADE_internal_data::ProSHADE_data::reportSymmetryResultsList
void reportSymmetryResultsList(ProSHADE_settings *settings)
This function takes prints the report for symmetry detection using multiple thresholds....
Definition: ProSHADE_data.cpp:3569
ProSHADE_internal_data::ProSHADE_data::zDimIndices
proshade_unsign zDimIndices
This is the size of the map cell z dimension in indices.
Definition: ProSHADE_data.hpp:67
ProSHADE_internal_data::ProSHADE_data::detectSymmetryFromAngleAxisSpace
void detectSymmetryFromAngleAxisSpace(ProSHADE_settings *settings)
This function runs the symmetry detection algorithms on this structure using the angle-axis space and...
Definition: ProSHADE_data.cpp:1954
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, proshade_signed messageShift=0)
General stdout message printing.
Definition: ProSHADE_messages.cpp:71
ProSHADE_internal_data::ProSHADE_data::yDimSize
proshade_single yDimSize
This is the size of the map cell y dimension in Angstroms.
Definition: ProSHADE_data.hpp:60
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:470
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:4124
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:4154
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:117