ProSHADE  0.7.6.2 (DEC 2021)
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 * > *axes, std::vector< std::vector< proshade_double > > *allCs, std::vector< proshade_double > *mapCOMShift)
 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, std::vector< std::vector< proshade_double > > *allCs, std::vector< proshade_double * > *axes, 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 673 of file ProSHADE_tasks.cpp.

674 {
675  //================================================ Are the any structures?
676  if ( settings->inputFiles.size () != 2 )
677  {
678  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." );
679  }
680 
681  //================================================ If centring is on, turn it off and report warning.
682  if ( settings->moveToCOM )
683  {
684  ProSHADE_internal_messages::printWarningMessage ( settings->verbose, "!!! ProSHADE WARNING !!! Map centring was requested, but makes no sense for overlay mode. Turning it off.", "WO00066" );
685  settings->moveToCOM = false;
686  }
687 
688  //================================================ Done
689  return ;
690 
691 }

◆ 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 580 of file ProSHADE_tasks.cpp.

581 {
582  //================================================ Are the any structures?
583  if ( settings->inputFiles.size () < 1 )
584  {
585  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." );
586  }
587 
588  //================================================ Is the axis tolerance set properly?
589  if ( settings->axisErrTolerance < 0.0 )
590  {
591  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." );
592  }
593 
594  //================================================ Done
595  return ;
596 
597 }

◆ 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 609 of file ProSHADE_tasks.cpp.

610 {
611  //================================================ Check the settings are complete and meaningful
612  checkOverlaySettings ( settings );
613 
614  //================================================ Initialise variables
615  proshade_double eulA, eulB, eulG, trsX, trsY, trsZ;
616 
617  //================================================ Create the data objects initially (this time without phase)
620 
621  //================================================ First, run without phase and find best rotation angles
622  settings->usePhase = false;
623  ProSHADE_internal_overlay::getOptimalRotation ( settings, staticStructure, movingStructure, &eulA, &eulB, &eulG );
624 
625  //================================================ Release memory
626  delete staticStructure;
627  delete movingStructure;
628 
629  //================================================ Create the data objects again (this time with phase)
630  staticStructure = new ProSHADE_internal_data::ProSHADE_data ( );
631  movingStructure = new ProSHADE_internal_data::ProSHADE_data ( );
632 
633  //================================================ Now, run with phase and find optimal translation
634  settings->usePhase = true;
635  settings->changeMapResolution = true;
636  ProSHADE_internal_overlay::getOptimalTranslation ( settings, staticStructure, movingStructure, &trsX, &trsY, &trsZ, eulA, eulB, eulG );
637 
638  //================================================ Compute the proper translations using the translation function output
639  ProSHADE_internal_misc::addToDoubleVector ( rotationCentre, movingStructure->originalPdbRotCenX );
640  ProSHADE_internal_misc::addToDoubleVector ( rotationCentre, movingStructure->originalPdbRotCenY );
641  ProSHADE_internal_misc::addToDoubleVector ( rotationCentre, movingStructure->originalPdbRotCenZ );
642  ProSHADE_internal_misc::addToDoubleVector ( finalTranslation, movingStructure->originalPdbTransX );
643  ProSHADE_internal_misc::addToDoubleVector ( finalTranslation, movingStructure->originalPdbTransY );
644  ProSHADE_internal_misc::addToDoubleVector ( finalTranslation, movingStructure->originalPdbTransZ );
645 
646  //================================================ Write out everything
647  movingStructure->writeOutOverlayFiles ( settings, eulA, eulB, eulG, rotationCentre, finalTranslation );
648 
649  //================================================ Save the rotation and rest of translations
650  ProSHADE_internal_misc::addToDoubleVector ( eulerAngles, eulA );
651  ProSHADE_internal_misc::addToDoubleVector ( eulerAngles, eulB );
652  ProSHADE_internal_misc::addToDoubleVector ( eulerAngles, eulG );
653 
654  //================================================ Report results to user
655  movingStructure->reportOverlayResults ( settings, rotationCentre, eulerAngles, finalTranslation );
656 
657  //================================================ Release memory
658  delete staticStructure;
659  delete movingStructure;
660 
661  //================================================ Done
662  return ;
663 
664 }

◆ 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,
std::vector< std::vector< proshade_double > > *  allCs,
std::vector< proshade_double * > *  axes,
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 366 of file ProSHADE_tasks.cpp.

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

◆ SymmetryDetectionTask()

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 
)

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]axesA pointer to a vector to which all the axes of the recommended symmetry (if any) will be saved.
[in]allCsA pointer to a vector to which all the detected cyclic symmetries will be saved into.
[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.

Definition at line 288 of file ProSHADE_tasks.cpp.

289 {
290  //================================================ Check the settings are complete and meaningful
291  checkSymmetrySettings ( settings );
292 
293  //================================================ Now, for each other structure
294  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( settings->inputFiles.size() ); iter++ )
295  {
296  //============================================ Create a data object
298 
299  //============================================ Read in the compared structure
300  symmetryStructure->readInStructure ( settings->inputFiles.at(iter), iter, settings );
301 
302  //============================================ Assume symmetry centre at the box centre, or find it out using Patterson map?
303  if ( settings->findSymCentre )
304  {
305  //======================================== Report progress
306  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 1, "Starting symmetry centre detection procedure.", settings->messageShift );
307 
308  //======================================== Make a local copy of settings (to avoid centre detection settings things for the symmetry detection which will follow)
309  ProSHADE_settings* rotCenSettings = new ProSHADE_settings ( settings );
310  rotCenSettings->messageShift = 1;
311 
312  //======================================== Run the detection
313  SymmetryCentreDetectionTask ( rotCenSettings, allCs, axes, iter );
314 
315  //======================================== Save the results
316  settings->centrePosition.at(0) = rotCenSettings->centrePosition.at(0);
317  settings->centrePosition.at(1) = rotCenSettings->centrePosition.at(1);
318  settings->centrePosition.at(2) = rotCenSettings->centrePosition.at(2);
319 
320  //======================================== Report progress
321  std::stringstream ss;
322  ss << "Detected symmetry centre at " << settings->centrePosition.at(0) << " ; " << settings->centrePosition.at(1) << " ; " << settings->centrePosition.at(2) << std::endl;
323  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, ss.str(), settings->messageShift );
324  }
325 
326  //============================================ Internal data processing (COM, norm, mask, extra space)
327  symmetryStructure->processInternalMap ( settings );
328 
329  //============================================ Map to sphere
330  symmetryStructure->mapToSpheres ( settings );
331 
332  //============================================ Get spherical harmonics
333  symmetryStructure->computeSphericalHarmonics ( settings );
334 
335  //============================================ Compute auto-rotation map
336  symmetryStructure->computeRotationFunction ( settings );
337 
338  //======================================== Detect point groups in the angle-axis space
339  symmetryStructure->detectSymmetryFromAngleAxisSpace ( settings, axes, allCs );
340 
341  //============================================ Report results
342  symmetryStructure->reportSymmetryResults ( settings );
343 
344  //============================================ Save internal map shift to run object,
345  ProSHADE_internal_misc::addToDoubleVector ( mapCOMShift, symmetryStructure->mapCOMProcessChangeX );
346  ProSHADE_internal_misc::addToDoubleVector ( mapCOMShift, symmetryStructure->mapCOMProcessChangeY );
347  ProSHADE_internal_misc::addToDoubleVector ( mapCOMShift, symmetryStructure->mapCOMProcessChangeZ );
348 
349  //============================================ Release memory
350  delete symmetryStructure;
351  }
352 
353  //================================================ Done
354  return ;
355 
356 }
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:1847
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:115
ProSHADE_settings::findSymCentre
bool findSymCentre
Should phase-less map be used to determine centre of symmetry?
Definition: ProSHADE_settings.hpp:126
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:116
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:1374
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:1316
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:1820
ProSHADE_internal_tasks::checkDistancesSettings
void checkDistancesSettings(ProSHADE_settings *settings)
The distances computation settings checks.
Definition: ProSHADE_tasks.cpp:258
ProSHADE_internal_symmetry::findReliableUnphasedSymmetries
std::vector< proshade_unsign > findReliableUnphasedSymmetries(std::vector< std::vector< proshade_double > > *allCs, 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:3674
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:110
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:3935
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:133
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:3791
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:908
ProSHADE_internal_data::ProSHADE_data
This class contains all inputed and derived data for a single structure.
Definition: ProSHADE_data.hpp:49
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_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)
This function finds the Centre of Mass for a map.
Definition: ProSHADE_mapManip.cpp:240
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:673
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: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:3227
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:3959
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:149
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:3447
ProSHADE_settings::messageShift
proshade_signed messageShift
This value allows shifting the messages to create more readable log for sub-processes.
Definition: ProSHADE_settings.hpp:150
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::getXDim
proshade_unsign getXDim(void)
This function allows access to the map size in indices along the X axis.
Definition: ProSHADE_data.cpp:3929
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:509
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:1803
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:1788
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:4512
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:892
ProSHADE_settings::reBoxMap
bool reBoxMap
This switch decides whether re-boxing is needed.
Definition: ProSHADE_settings.hpp:92
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:103
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:3939
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:4478
ProSHADE_internal_tasks::SymmetryCentreDetectionTask
void SymmetryCentreDetectionTask(ProSHADE_settings *settings, std::vector< std::vector< proshade_double > > *allCs, std::vector< proshade_double * > *axes, proshade_unsign strIndex=0)
The task for finding the structure centre based on phase-less symmetry detection.
Definition: ProSHADE_tasks.cpp:366
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:1695
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:3949
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:142
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::getMapValue
proshade_double getMapValue(proshade_unsign pos)
This function returns the internal map representation value of a particular array position.
Definition: ProSHADE_data.cpp:3601
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:3831
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:3999
ProSHADE_settings::moveToCOM
bool moveToCOM
Logical value stating whether the structure should be moved to have its Centre Of Mass (COM) in the m...
Definition: ProSHADE_settings.hpp:99
ProSHADE_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:3989
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:3969
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:3479
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:128
ProSHADE_internal_tasks::checkSymmetrySettings
void checkSymmetrySettings(ProSHADE_settings *settings)
The symmetry computation settings checks.
Definition: ProSHADE_tasks.cpp:580
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, 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:1907
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_data::ProSHADE_data::getAllGroupElements
std::vector< std::vector< proshade_double > > getAllGroupElements(ProSHADE_settings *settings, 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:3151
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:372
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:3979
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:4009
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:113