ProSHADE  0.7.6.0 (JUL 2021)
Protein Shape Detection
ProSHADE_distances.cpp
Go to the documentation of this file.
1 
23 //==================================================== ProSHADE
24 #include "ProSHADE_distances.hpp"
25 
32 {
33  //================================================ Allocate the required memory
34  this->rrpMatrices = new proshade_double** [this->maxShellBand];
35  ProSHADE_internal_misc::checkMemoryAllocation ( this->rrpMatrices, __FILE__, __LINE__, __func__ );
36 
37  for ( proshade_unsign bwIt = 0; bwIt < this->maxShellBand; bwIt++ )
38  {
39  //============================================ For rach sphere
40  this->rrpMatrices[bwIt] = new proshade_double* [this->noSpheres];
41  ProSHADE_internal_misc::checkMemoryAllocation ( this->rrpMatrices[bwIt], __FILE__, __LINE__, __func__ );
42 
43  for ( proshade_unsign shIt = 0; shIt < this->noSpheres; shIt++ )
44  {
45  this->rrpMatrices[bwIt][shIt] = new double [this->noSpheres];
46  ProSHADE_internal_misc::checkMemoryAllocation ( this->rrpMatrices[bwIt][shIt], __FILE__, __LINE__, __func__ );
47  }
48  }
49 }
50 
60 {
61  //================================================ Report progress
62  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, "Computing RRP matrices for structure " + this->fileName );
63 
64  //================================================ Allocate the memory
65  this->allocateRRPMemory ( );
66 
67  //================================================ Start computation: For each band (l)
68  proshade_double descValR = 0.0;
69  proshade_unsign arrPos1, arrPos2;
70  for ( proshade_unsign band = 0; band < this->maxShellBand; band++ )
71  {
72  //============================================ For each unique shell couple
73  for ( proshade_unsign shell1 = 0; shell1 < this->noSpheres; shell1++ )
74  {
75  //======================================== Does the band exist for this shell1?
76  if ( !ProSHADE_internal_distances::isBandWithinShell ( band, shell1, this->spheres ) )
77  {
78  for ( proshade_unsign shell2 = 0; shell2 < this->noSpheres; shell2++ )
79  {
80  this->rrpMatrices[band][shell1][shell2] = 0.0;
81  this->rrpMatrices[band][shell2][shell1] = 0.0;
82  }
83  continue;
84  }
85 
86  for ( proshade_unsign shell2 = 0; shell2 < this->noSpheres; shell2++ )
87  {
88  //==================================== Compute each values only once
89  if ( shell1 > shell2 ) { continue; }
90 
91  //==================================== Check if band exists for this shell2?
92  if ( !ProSHADE_internal_distances::isBandWithinShell ( band, shell2, this->spheres ) )
93  {
94  this->rrpMatrices[band][shell1][shell2] = 0.0;
95  this->rrpMatrices[band][shell2][shell1] = 0.0;
96  continue;
97  }
98 
99  //==================================== Initialise
100  descValR = 0.0;
101 
102  //==================================== Sum over order (m)
103  for ( proshade_unsign order = 0; order < static_cast< proshade_unsign > ( ( 2 * band ) + 1 ); order++ )
104  {
105  arrPos1 = static_cast< proshade_unsign > ( seanindex ( static_cast< int > ( order ) - static_cast<int > ( band ),
106  static_cast< int > ( band ), static_cast< int > ( this->spheres[shell1]->getLocalBandwidth() ) ) );
107  arrPos2 = static_cast< proshade_unsign > ( seanindex ( static_cast< int > ( order ) - static_cast< int > ( band ),
108  static_cast< int > ( band ), static_cast< int > ( this->spheres[shell2]->getLocalBandwidth() ) ) );
109  descValR += ProSHADE_internal_maths::complexMultiplicationConjugRealOnly ( &this->sphericalHarmonics[shell1][arrPos1][0],
110  &this->sphericalHarmonics[shell1][arrPos1][1],
111  &this->sphericalHarmonics[shell2][arrPos2][0],
112  &this->sphericalHarmonics[shell2][arrPos2][1] );
113  }
114 
115  //==================================== Save the matrices
116  this->rrpMatrices[band][shell1][shell2] = descValR;
117  this->rrpMatrices[band][shell2][shell1] = descValR;
118  }
119  }
120  }
121 
122  //================================================ Report progress
123  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 3, "RRP matrices successfully computed." );
124 
125  //================================================ Done
126  return ;
127 
128 }
129 
139 bool ProSHADE_internal_distances::isBandWithinShell ( proshade_unsign bandInQuestion, proshade_unsign shellInQuestion, ProSHADE_internal_spheres::ProSHADE_sphere** spheres )
140 {
141  if ( bandInQuestion < spheres[shellInQuestion]->getLocalBandwidth() )
142  {
143  return ( true );
144  }
145  else
146  {
147  return ( false );
148  }
149 }
150 
162 {
163  //================================================ Report starting the task
164  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 1, "Starting energy levels distance computation." );
165 
166  //================================================ Initialise local variables
167  proshade_double ret = 0.0;
168  std::vector<proshade_double> bandDists;
169 
170  //================================================ Sanity check
171  if ( !settings->computeEnergyLevelsDesc )
172  {
173  throw ProSHADE_exception ( "Attempted computing energy levels descriptors when it was not required.", "ED00017", __FILE__, __LINE__, __func__, "Attempted to pre-compute the RRP matrices, when the user\n : has specifically stated that these should not be computed.\n : Unless you manipulated the code, this error should never\n : occur; if you see this, I made a large blunder. Please let\n : me know!" );
174  }
175 
176  //================================================ Get the RRP matrices for both objects
177  obj1->computeRRPMatrices ( settings );
178  obj2->computeRRPMatrices ( settings );
179 
180  //================================================ Find the minimium comparable shells and bands
181  proshade_unsign minCommonShells = std::min ( obj1->getMaxSpheres(), obj2->getMaxSpheres() );
182  proshade_unsign minCommonBands = std::min ( obj1->getMaxBand(), obj2->getMaxBand() );
183 
184  //================================================ Get the Pearson's coefficients for each common band
185  computeRRPPearsonCoefficients ( obj1, obj2, settings, minCommonBands, minCommonShells, &bandDists );
186 
187  //================================================ Get distance (by averaging Patterson's coefficients)
188  ret = static_cast<proshade_double> ( std::accumulate ( bandDists.begin(), bandDists.end(), 0.0 ) ) /
189  static_cast<proshade_double> ( bandDists.size() );
190 
191  //================================================ Report completion
192  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, "Energy levels distance computation complete." );
193 
194  //================================================ Done
195  return ( ret );
196 
197 }
198 
211 void ProSHADE_internal_distances::computeRRPPearsonCoefficients ( ProSHADE_internal_data::ProSHADE_data* obj1, ProSHADE_internal_data::ProSHADE_data* obj2, ProSHADE_settings* settings, proshade_unsign minCommonBands, proshade_unsign minCommonShells, std::vector<proshade_double>* bandDists )
212 {
213  //================================================ Report completion
214  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, "Correlating RRP matrices." );
215 
216  //================================================ Initialise local variables
217  proshade_double *str1Vals = new proshade_double[minCommonShells * minCommonShells];
218  proshade_double *str2Vals = new proshade_double[minCommonShells * minCommonShells];
219  ProSHADE_internal_misc::checkMemoryAllocation ( str1Vals, __FILE__, __LINE__, __func__ );
220  ProSHADE_internal_misc::checkMemoryAllocation ( str2Vals, __FILE__, __LINE__, __func__ );
221  proshade_unsign arrIter = 0;
222 
223  //================================================ Start computation: For each band (l)
224  for ( proshade_unsign band = 0; band < minCommonBands; band++ )
225  {
226  //============================================ Reset local counter
227  arrIter = 0;
228 
229  //============================================ For each shell pair
230  for ( proshade_unsign shell1 = 0; shell1 < minCommonShells; shell1++ )
231  {
232  //======================================== Check if band exists (progressive only)
233  if ( settings->progressiveSphereMapping ) { if ( !obj1->shellBandExists( shell1, band ) || !obj2->shellBandExists( shell1, band ) ) { continue; } }
234 
235  for ( proshade_unsign shell2 = 0; shell2 < minCommonShells; shell2++ )
236  {
237  //============================ Check the other shell as well
238  if ( !obj1->shellBandExists( shell2, band ) || !obj2->shellBandExists( shell2, band ) ) { continue; }
239 
240  //==================================== Set values between which the Person's correlation coefficient should be computed
241  str1Vals[arrIter] = obj1->getRRPValue ( band, shell1, shell2 ) *
242  pow ( static_cast<proshade_double> ( shell1 ), settings->enLevMatrixPowerWeight ) *
243  pow ( static_cast<proshade_double> ( shell2 ), settings->enLevMatrixPowerWeight );
244  str2Vals[arrIter] = obj2->getRRPValue ( band, shell1, shell2 ) *
245  pow ( static_cast<proshade_double> ( shell1 ), settings->enLevMatrixPowerWeight ) *
246  pow ( static_cast<proshade_double> ( shell2 ), settings->enLevMatrixPowerWeight );
247 
248  arrIter += 1;
249  }
250  }
251 
252  //============================================ Get Pearson's Correlation Coefficient
253  ProSHADE_internal_misc::addToDoubleVector ( bandDists, ProSHADE_internal_maths::pearsonCorrCoeff ( str1Vals, str2Vals, arrIter ) );
254  }
255 
256  //================================================ Clean up
257  delete[] str1Vals;
258  delete[] str2Vals;
259 
260  //================================================ Report completion
261  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 3, "RRP matrices correlation computed." );
262 
263  //================================================ Done
264  return ;
265 }
266 
276 {
277  //================================================ Save the maximum band to the object
278  this->maxCompBand = band;
279 
280  //================================================ Allocate the required memory
281  this->eMatrices = new proshade_complex** [this->maxCompBand];
282  ProSHADE_internal_misc::checkMemoryAllocation ( this->eMatrices, __FILE__, __LINE__, __func__ );
283 
284  for ( proshade_unsign bandIter = 0; bandIter < this->maxCompBand; bandIter++ )
285  {
286  //============================================ Allocate the data structure
287  this->eMatrices[bandIter] = new proshade_complex* [static_cast<proshade_unsign> ( ( bandIter * 2 ) + 1 )];
288  ProSHADE_internal_misc::checkMemoryAllocation ( this->eMatrices[bandIter], __FILE__, __LINE__, __func__ );
289 
290  for ( proshade_unsign band2Iter = 0; band2Iter < static_cast<proshade_unsign> ( ( bandIter * 2 ) + 1 ); band2Iter++ )
291  {
292  this->eMatrices[bandIter][band2Iter] = new proshade_complex [static_cast<proshade_unsign> ( ( bandIter * 2 ) + 1 )];
293  ProSHADE_internal_misc::checkMemoryAllocation ( this->eMatrices[bandIter][band2Iter], __FILE__, __LINE__, __func__ );
294  }
295  }
296 
297  //================================================ Done
298  return ;
299 
300 }
301 
312 void ProSHADE_internal_distances::allocateTrSigmaWorkspace ( proshade_unsign minSpheres, proshade_unsign intOrder, proshade_double*& obj1Vals, proshade_double*& obj2Vals, proshade_double*& GLabscissas, proshade_double*& GLweights, proshade_complex*& radiiVals )
313 {
314  //================================================ Allocate the memory
315  obj1Vals = new proshade_double [minSpheres];
316  obj2Vals = new proshade_double [minSpheres];
317  radiiVals = new proshade_complex[minSpheres];
318  GLabscissas = new proshade_double [intOrder];
319  GLweights = new proshade_double [intOrder];
320 
321  //================================================ Check the memory allocation
322  ProSHADE_internal_misc::checkMemoryAllocation ( obj1Vals, __FILE__, __LINE__, __func__ );
323  ProSHADE_internal_misc::checkMemoryAllocation ( obj2Vals, __FILE__, __LINE__, __func__ );
324  ProSHADE_internal_misc::checkMemoryAllocation ( radiiVals, __FILE__, __LINE__, __func__ );
325  ProSHADE_internal_misc::checkMemoryAllocation ( GLabscissas, __FILE__, __LINE__, __func__ );
326  ProSHADE_internal_misc::checkMemoryAllocation ( GLweights, __FILE__, __LINE__, __func__ );
327 
328  //================================================ Done
329  return ;
330 
331 }
332 
341 void ProSHADE_internal_distances::computeSphericalHarmonicsMagnitude ( ProSHADE_internal_data::ProSHADE_data* obj, proshade_unsign band, proshade_unsign order, proshade_unsign radius, proshade_double* result )
342 {
343  //================================================ Find the magnitude
345  obj->getImagSphHarmValue ( band, order, radius ),
346  obj->getRealSphHarmValue ( band, order, radius ),
347  obj->getImagSphHarmValue ( band, order, radius ) );
348 
349  //================================================ Weight by radius^2 for the integration that will follow
350  *result *= pow ( static_cast<proshade_double> ( obj->getAnySphereRadius( radius ) ), 2.0 );
351 
352  //================================================ Done
353  return ;
354 
355 }
356 
370 void ProSHADE_internal_distances::computeEMatricesForLM ( ProSHADE_internal_data::ProSHADE_data* obj1, ProSHADE_internal_data::ProSHADE_data* obj2, proshade_unsign bandIter, proshade_unsign orderIter, proshade_complex* radiiVals, proshade_unsign integOrder, proshade_double* abscissas, proshade_double* weights, proshade_double integRange, proshade_double sphereDist )
371 {
372  //================================================ Initialise local variables
373  proshade_unsign objCombValsIter = 0;
374  proshade_double hlpReal, hlpImag;
375  proshade_complex arrVal;
376 
377  //================================================ For each combination of m and m' for E matrices
378  for ( proshade_unsign order2Iter = 0; order2Iter < ( ( bandIter * 2 ) + 1 ); order2Iter++ )
379  {
380  //============================================ Reset loop
381  objCombValsIter = 0;
382 
383  //============================================ Find the c*conj(c) values for different radii
384  for ( proshade_unsign radiusIter = 0; radiusIter < std::min( obj1->getMaxSpheres(), obj2->getMaxSpheres() ); radiusIter++ )
385  {
386 
387  //======================================== Get only values where the shell has the band
388  if ( std::min ( obj1->getShellBandwidth ( radiusIter ), obj2->getShellBandwidth ( radiusIter ) ) <= bandIter ) { continue; }
389 
390  //======================================== Multiply coeffs
391  ProSHADE_internal_maths::complexMultiplicationConjug ( obj1->getRealSphHarmValue ( bandIter, orderIter, radiusIter ),
392  obj1->getImagSphHarmValue ( bandIter, orderIter, radiusIter ),
393  obj2->getRealSphHarmValue ( bandIter, order2Iter, radiusIter ),
394  obj2->getImagSphHarmValue ( bandIter, order2Iter, radiusIter ),
395  &hlpReal, &hlpImag );
396 
397  //======================================== Apply r^2 integral weight
398  radiiVals[objCombValsIter][0] = hlpReal * pow ( ( static_cast<proshade_double> ( obj1->getAnySphereRadius( radiusIter ) ) ), 2.0 );
399  radiiVals[objCombValsIter][1] = hlpImag * pow ( ( static_cast<proshade_double> ( obj1->getAnySphereRadius( radiusIter ) ) ), 2.0 );
400 
401  objCombValsIter += 1;
402  }
403 
404  //============================================ Integrate over all radii using n-point Gauss-Legendre integration
405  ProSHADE_internal_maths::gaussLegendreIntegration ( radiiVals, objCombValsIter, integOrder, abscissas, weights, integRange, sphereDist, &hlpReal, &hlpImag );
406 
407  //============================================ Save the result into E matrices
408  arrVal[0] = hlpReal;
409  arrVal[1] = hlpImag;
410  obj2->setEMatrixValue ( bandIter, orderIter, order2Iter, arrVal );
411  }
412 
413  //================================================ Done
414  return ;
415 
416 }
417 
433 proshade_double ProSHADE_internal_distances::computeWeightsForEMatricesForLM ( ProSHADE_internal_data::ProSHADE_data* obj1, ProSHADE_internal_data::ProSHADE_data* obj2, proshade_unsign bandIter, proshade_unsign orderIter, proshade_double* obj1Vals, proshade_double* obj2Vals, proshade_unsign integOrder, proshade_double* abscissas, proshade_double* weights, proshade_single sphereDist )
434 {
435  //================================================ Initialise local values
436  proshade_unsign obj1ValsIter = 0;
437  proshade_unsign obj2ValsIter = 0;
438 
439  //================================================ Set sphere counters
440  proshade_unsign minSphere = std::min( obj1->getMaxSpheres(), obj2->getMaxSpheres() );
441  proshade_unsign maxSphere = 0;
442 
443  //================================================ For each radius, deal with weights
444  for ( proshade_unsign radiusIter = 0; radiusIter < std::min( obj1->getMaxSpheres(), obj2->getMaxSpheres() ); radiusIter++ )
445  {
446  //============================================ Get only values where the shell has the band
447  if ( std::min ( obj1->getShellBandwidth ( radiusIter ), obj2->getShellBandwidth ( radiusIter ) ) <= bandIter ) { continue; }
448  minSphere = std::min ( radiusIter, minSphere );
449  maxSphere = std::max ( radiusIter, maxSphere );
450 
451  //============================================ Get the magnitudes for weighting
452  computeSphericalHarmonicsMagnitude ( obj1, bandIter, orderIter, radiusIter, &(obj1Vals[obj1ValsIter]) );
453  computeSphericalHarmonicsMagnitude ( obj2, bandIter, orderIter, radiusIter, &(obj2Vals[obj2ValsIter]) );
454  obj1ValsIter += 1;
455  obj2ValsIter += 1;
456  }
457 
458  //================================================ Integrate weights
459  proshade_single minSphereRad = obj1->getSpherePosValue ( minSphere ) - ( sphereDist * 0.5f );
460  proshade_single maxSphereRad = obj1->getSpherePosValue ( maxSphere ) + ( sphereDist * 0.5f );
461 
462  obj1->setIntegrationWeightCumul ( ProSHADE_internal_maths::gaussLegendreIntegrationReal ( obj1Vals, obj1ValsIter, integOrder, abscissas, weights, static_cast< proshade_double > ( maxSphereRad - minSphereRad ), static_cast< proshade_double > ( sphereDist ) ) );
463  obj2->setIntegrationWeightCumul ( ProSHADE_internal_maths::gaussLegendreIntegrationReal ( obj2Vals, obj2ValsIter, integOrder, abscissas, weights, static_cast< proshade_double > ( maxSphereRad - minSphereRad ), static_cast< proshade_double > ( sphereDist ) ) );
464 
465  //================================================ Done
466  return ( static_cast< proshade_double > ( maxSphereRad - minSphereRad ) );
467 
468 }
469 
478 void ProSHADE_internal_distances::releaseTrSigmaWorkspace ( proshade_double*& obj1Vals, proshade_double*& obj2Vals, proshade_double*& GLabscissas, proshade_double*& GLweights, proshade_complex*& radiiVals )
479 {
480  //================================================ Release memory
481  delete[] obj1Vals;
482  delete[] obj2Vals;
483  delete[] radiiVals;
484  delete[] GLabscissas;
485  delete[] GLweights;
486 
487  //================================================ Set to NULL
488  obj1Vals = nullptr;
489  obj2Vals = nullptr;
490  radiiVals = nullptr;
491  GLabscissas = nullptr;
492  GLweights = nullptr;
493 
494  //================================================ Done
495  return ;
496 
497 }
498 
511 {
512  //================================================ Report progress
513  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, "Starting computation of E matrices." );
514 
515  //================================================ Allocatre memory for E matrices in the second object (first may be compared to more structures and therefore its data would be written over)
516  obj2->allocateEMatrices ( std::min ( obj1->getMaxBand(), obj2->getMaxBand() ) );
517 
518  //================================================ Initialise local variables
519  proshade_double *obj1Vals, *obj2Vals, *GLAbscissas, *GLWeights;
520  proshade_complex* radiiVals;
521  proshade_double integRange;
522 
523  //================================================ Allocate workspace memory
524  allocateTrSigmaWorkspace ( std::min( obj1->getMaxSpheres(), obj2->getMaxSpheres() ), settings->integOrder, obj1Vals, obj2Vals, GLAbscissas, GLWeights, radiiVals);
525 
526  //================================================ Initialise abscissas and weights for integration
527  ProSHADE_internal_maths::getLegendreAbscAndWeights ( settings->integOrder, GLAbscissas, GLWeights, settings->taylorSeriesCap );
528 
529  //================================================ For each band (l), compute the E matrix integrals
530  for ( proshade_unsign bandIter = 0; bandIter < std::min ( obj1->getMaxBand(), obj2->getMaxBand() ); bandIter++ )
531  {
532  //============================================ For each order (m)
533  for ( proshade_unsign orderIter = 0; orderIter < ( ( bandIter * 2 ) + 1 ); orderIter++ )
534  {
535  //======================================== Get weights for the required band(l) and order (m)
536  integRange = computeWeightsForEMatricesForLM ( obj1, obj2, bandIter, orderIter, obj1Vals, obj2Vals, settings->integOrder, GLAbscissas, GLWeights, settings->maxSphereDists );
537 
538  //======================================== Compute E matrices value for given band (l) and order(m)
539  computeEMatricesForLM ( obj1, obj2, bandIter, orderIter, radiiVals, settings->integOrder, GLAbscissas, GLWeights, integRange, static_cast< proshade_double > ( settings->maxSphereDists ) );
540  }
541 
542  //============================================ Report progress
543  if ( settings->verbose > 3 )
544  {
545  std::stringstream hlpSS;
546  hlpSS << "E matrices computed for band " << bandIter;
547  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 4, hlpSS.str() );
548  }
549  }
550 
551  //================================================ Release the workspace memory
552  releaseTrSigmaWorkspace ( obj1Vals, obj2Vals, GLAbscissas, GLWeights, radiiVals );
553 
554  //================================================ Report progress
555  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 3, "E matrices computed." );
556 
557  //================================================ Done
558  return ;
559 
560 }
561 
573 {
574  //================================================ Report progress
575  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 3, "Starting E matrices normalisation." );
576 
577  //================================================ Normalise by the Pearson's c.c. like formula
578  proshade_double eMatNormFactor = std::sqrt ( obj1->getIntegrationWeight() * obj2->getIntegrationWeight() );
579 
580  //================================================ If this is self-correlation (i.e. obj1 == obj2), then divide normalisation factor by 2 as the weight was applied cumulatively!
581  if ( settings->task == Symmetry ) { eMatNormFactor /= 2.0; }
582 
583  for ( proshade_unsign bandIter = 0; bandIter < std::min ( obj1->getMaxBand(), obj2->getMaxBand() ); bandIter++ )
584  {
585  //============================================ For each combination of m and m' for E matrices
586  for ( proshade_unsign orderIter = 0; orderIter < ( ( bandIter * 2 ) + 1 ); orderIter++ )
587  {
588  for ( proshade_unsign order2Iter = 0; order2Iter < ( ( bandIter * 2 ) + 1 ); order2Iter++ )
589  {
590  obj2->normaliseEMatrixValue ( bandIter, orderIter, order2Iter, eMatNormFactor );
591  }
592  }
593  }
594 
595  //================================================ Report progress
596  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 4, "E matrices normalised." );
597 
598  //================================================ Done
599  return ;
600 
601 }
602 
617 {
618  //================================================ Report starting the task
619  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 1, "Starting trace sigma distance computation." );
620 
621  //================================================ Initialise return variable
622  proshade_double ret = 0.0;
623 
624  //================================================ Sanity check
625  if ( !settings->computeTraceSigmaDesc )
626  {
627  throw ProSHADE_exception ( "Attempted computing trace sigma descriptors when it was\n : not required.", "ED00018", __FILE__, __LINE__, __func__, "Attempted to pre-compute the E matrices, when the user\n : has specifically stated that these should not be computed.\n : Unless you manipulated the code, this error should never\n : occur; if you see this, I made a large blunder. Please let\n : me know!" );
628  }
629 
630  //================================================ Empty the cumulative weights back to 0.0 for each structure
631  obj1->setIntegrationWeight ( 0.0 );
632  obj1->setIntegrationWeight ( 0.0 );
633 
634  //================================================ Compute un-weighted E matrices and their weights
635  computeEMatrices ( obj1, obj2, settings );
636 
637  //================================================ Normalise E matrices by the magnitudes
638  normaliseEMatrices ( obj1, obj2, settings );
639 
640  //================================================ Allocate the required memory
641  double* singularValues = new double[( ( std::min ( obj1->getMaxBand(), obj2->getMaxBand() ) * 2 ) + 1 )];
642  ProSHADE_internal_misc::checkMemoryAllocation ( singularValues, __FILE__, __LINE__, __func__ );
643 
644  //================================================ Compute the distance
645  for ( proshade_unsign lIter = 0; lIter < std::min ( obj1->getMaxBand(), obj2->getMaxBand() ); lIter++ )
646  {
647  //============================================ Find the complex matrix SVD singular values
648  ProSHADE_internal_maths::complexMatrixSVDSigmasOnly ( obj2->getEMatrixByBand ( lIter ), static_cast<int> ( ( lIter * 2 ) + 1 ), singularValues );
649 
650  //============================================ Now sum the trace
651  for ( proshade_unsign iter = 0; iter < ( ( lIter * 2 ) + 1 ); iter++ )
652  {
653  ret += singularValues[iter];
654  }
655  }
656 
657  //================================================ Report completion
658  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 3, "E matrices decomposed to singular values." );
659 
660  //================================================ Release the memory
661  delete[] singularValues;
662 
663  //================================================ Report completion
664  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, "Trace sigma distance computation complete." );
665 
666  //================================================ Done
667  return ( ret );
668 
669 }
670 
676 {
677  //================================================ Allocate the memory
678  this->so3Coeffs = new fftw_complex [static_cast<proshade_unsign>( ( 4 * pow( static_cast<proshade_double> ( band ), 3.0 ) - static_cast<proshade_double> ( band ) ) / 3.0 )];
679  this->so3CoeffsInverse = new fftw_complex [static_cast<proshade_unsign>( pow( static_cast<proshade_double> ( band ) * 2.0, 3.0 ) )];
680 
681  //================================================ Check memory allocation
682  ProSHADE_internal_misc::checkMemoryAllocation ( this->so3Coeffs, __FILE__, __LINE__, __func__ );
683  ProSHADE_internal_misc::checkMemoryAllocation ( this->so3CoeffsInverse, __FILE__, __LINE__, __func__ );
684 
685  //================================================ Done
686  return ;
687 
688 }
689 
702 {
703  //================================================ Report progress
704  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, "Converting E matrices to SO(3) coefficients." );
705 
706  //================================================ Allocate memory for the coefficients
707  obj2->allocateSO3CoeffsSpace ( std::min ( obj1->getMaxBand(), obj2->getMaxBand() ) );
708 
709  //================================================ Initialise local variables
710  proshade_double wigNorm, hlpValReal, hlpValImag;
711  proshade_double signValue = 1.0;
712  proshade_unsign indexO;
713  proshade_complex hlpVal;
714 
715  //================================================ For each band (l)
716  for ( proshade_signed bandIter = 0; bandIter < static_cast<proshade_signed> ( std::min ( obj1->getMaxBand(), obj2->getMaxBand() ) ); bandIter++ )
717  {
718  //============================================ Get wigner normalisation factor
719  wigNorm = 2.0 * M_PI * sqrt ( 2.0 / (2.0 * static_cast< proshade_double > ( bandIter ) + 1.0 ) ) ;
720 
721  //============================================ For each order (m)
722  for ( proshade_signed orderIter = 0; orderIter < ( ( bandIter * 2 ) + 1 ); orderIter++ )
723  {
724  //======================================== Set the sign
725  if ( ( orderIter - bandIter + bandIter ) % 2 ) { signValue = -1.0 ; }
726  else { signValue = 1.0 ; }
727 
728  //======================================== For each order2 (m')
729  for ( proshade_signed order2Iter = 0; order2Iter < ( ( bandIter * 2 ) + 1 ); order2Iter++ )
730  {
731  //==================================== Find output index
732  indexO = static_cast< proshade_unsign > ( so3CoefLoc ( static_cast< int > ( orderIter - bandIter ), static_cast< int > ( order2Iter - bandIter ), static_cast< int > ( bandIter ), static_cast< int > ( std::min ( obj1->getMaxBand(), obj2->getMaxBand() ) ) ) );
733 
734  //==================================== Compute and save the SO(3) coefficients
735  obj2->getEMatrixValue ( static_cast< proshade_unsign > ( bandIter ), static_cast< proshade_unsign > ( orderIter ), static_cast< proshade_unsign > ( order2Iter ), &hlpValReal, &hlpValImag );
736  hlpVal[0] = hlpValReal * wigNorm * signValue;
737  hlpVal[1] = hlpValImag * wigNorm * signValue;
738  obj2->setSO3CoeffValue ( indexO, hlpVal );
739 
740  //==================================== Switch the sign value
741  signValue *= -1.0;
742  }
743  }
744  }
745 
746  //================================================ Report progress
747  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 3, "SO(3) coefficients obtained." );
748 
749  //================================================ Done
750  return ;
751 
752 }
753 
761 void ProSHADE_internal_distances::allocateInvSOFTWorkspaces ( proshade_complex*& work1, proshade_complex*& work2, proshade_double*& work3, proshade_unsign band )
762 {
763  //================================================ Allocate memory
764  work1 = new proshade_complex[8 * static_cast<proshade_unsign> ( pow( static_cast<double> ( band ), 3.0 ) )];
765  work2 = new proshade_complex[14 * static_cast<proshade_unsign> ( pow( static_cast<double> ( band ), 2.0 ) ) + (48 * band)];
766  work3 = new proshade_double [2 * static_cast<proshade_unsign> ( pow( static_cast<double> ( band ), 2.0 ) ) + (24 * band)];
767 
768  //================================================ Check the memory allocation
769  ProSHADE_internal_misc::checkMemoryAllocation ( work1, __FILE__, __LINE__, __func__ );
770  ProSHADE_internal_misc::checkMemoryAllocation ( work2, __FILE__, __LINE__, __func__ );
771  ProSHADE_internal_misc::checkMemoryAllocation ( work3, __FILE__, __LINE__, __func__ );
772 
773  //================================================ Done
774  return ;
775 
776 }
777 
785 void ProSHADE_internal_distances::prepareInvSOFTPlan ( fftw_plan* inverseSO3, int band, fftw_complex* work1, proshade_complex* invCoeffs )
786 {
787  //================================================ Prepare the plan describing variables
788  int howmany = 4 * band * band;
789  int idist = 2 * band;
790  int odist = 2 * band;
791  int rank = 2;
792 
793  int inembed[2], onembed[2];
794  inembed[0] = 2 * band;
795  inembed[1] = 4 * band * band;
796  onembed[0] = 2 * band;
797  onembed[1] = 4 * band * band;
798 
799  int istride = 1;
800  int ostride = 1;
801 
802  int na[2];
803  na[0] = 1;
804  na[1] = 2 * band;
805 
806  //================================================ Create the plan
807  *inverseSO3 = fftw_plan_many_dft ( rank,
808  na,
809  howmany,
810  work1,
811  inembed,
812  istride,
813  idist,
814  invCoeffs,
815  onembed,
816  ostride,
817  odist,
818  FFTW_FORWARD,
819  FFTW_ESTIMATE );
820 
821  //================================================ Done
822  return ;
823 
824 }
825 
832 void ProSHADE_internal_distances::releaseInvSOFTMemory ( proshade_complex*& work1, proshade_complex*& work2, proshade_double*& work3 )
833 {
834  //================================================ Release memory
835  delete[] work1;
836  delete[] work2;
837  delete[] work3;
838 
839  //================================================ Done
840  return ;
841 }
842 
855 {
856  //================================================ Report progress
857  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, "Computing inverse SO(3) Fourier transform." );
858 
859  //================================================ Initialise local variables
860  proshade_complex *workspace1, *workspace2;
861  proshade_double *workspace3;
862  fftw_plan inverseSO3;
863 
864  //================================================ Allocate memory for the workspaces
865  allocateInvSOFTWorkspaces ( workspace1, workspace2, workspace3, std::min ( obj1->getMaxBand(), obj2->getMaxBand() ) );
866 
867  //================================================ Prepare the FFTW plan
868  prepareInvSOFTPlan ( &inverseSO3, static_cast< int > ( std::min ( obj1->getMaxBand(), obj2->getMaxBand() ) ), workspace1, obj2->getInvSO3Coeffs ( ) );
869 
870  //================================================ Compute the transform
871  Inverse_SO3_Naive_fftw ( static_cast< int > ( std::min ( obj1->getMaxBand(), obj2->getMaxBand() ) ),
872  obj2->getSO3Coeffs ( ),
873  obj2->getInvSO3Coeffs ( ),
874  workspace1,
875  workspace2,
876  workspace3,
877  &inverseSO3,
878  0 );
879 
880  //================================================ Release memory
881  releaseInvSOFTMemory ( workspace1, workspace2, workspace3 );
882  fftw_destroy_plan ( inverseSO3 );
883 
884  //================================================ Report progress
885  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 3, "Inverse SO(3) Fourier transform computed." );
886 
887  //================================================ Done
888  return ;
889 
890 }
891 
907 {
908  //================================================ Report starting the task
909  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 1, "Starting rotation function distance computation." );
910 
911  //================================================ Initialise return variable
912  proshade_double ret = 0.0;
913  proshade_double eulA, eulB, eulG, EMatR, EMatI, WigDR, WigDI;
914 
915  //================================================ Sanity check
916  if ( !settings->computeRotationFuncDesc )
917  {
918  throw ProSHADE_exception ( "Attempted computing rotation function descriptors when it\n : was not required.", "ED00023", __FILE__, __LINE__, __func__, "Attempted to compute the SO(3) transform and the rotation \n : function descriptor when the user did not request this. \n : Unless you manipulated the code, this error should never \n : occur; if you see this, I made a large blunder. \n : Please let me know!" );
919  }
920 
921  //================================================ Compute weighted E matrices if not already present
922  if ( !settings->computeTraceSigmaDesc )
923  {
924  computeEMatrices ( obj1, obj2, settings );
925  normaliseEMatrices ( obj1, obj2, settings );
926  }
927 
928  //================================================ Generate SO(3) coefficients
929  generateSO3CoeffsFromEMatrices ( obj1, obj2, settings );
930 
931  //================================================ Compute the inverse SO(3) Fourier Transform (SOFT) on the newly computed coefficients
932  computeInverseSOFTTransform ( obj1, obj2, settings );
933 
934  //================================================ Get inverse SO(3) map top peak Euler angle values
936  std::min ( obj1->getMaxBand(), obj2->getMaxBand() ) * 2,
937  &eulA, &eulB, &eulG, settings );
938 
939  //================================================ Compute the Wigner D matrices for the Euler angles
940  ProSHADE_internal_wigner::computeWignerMatricesForRotation ( settings, obj2, eulA, eulB, eulG );
941 
942  //================================================ Compute the distance
943  for ( proshade_unsign bandIter = 0; bandIter < obj2->getComparisonBand(); bandIter++ )
944  {
945  //============================================ For each order1
946  for ( proshade_unsign order1 = 0; order1 < ( ( bandIter * 2 ) + 1 ); order1++ )
947  {
948  //======================================== For each order2
949  for ( proshade_unsign order2 = 0; order2 < ( ( bandIter * 2 ) + 1 ); order2++ )
950  {
951  //==================================== Multiply D_{l} * E_{l} and get sum over l of traces (i.e. just sum all together)
952  obj2->getEMatrixValue ( bandIter, order1, order2, &EMatR, &EMatI );
953  obj2->getWignerMatrixValue ( bandIter, order2, order1, &WigDR, &WigDI );
954  ret += ProSHADE_internal_maths::complexMultiplicationRealOnly ( &WigDR, &WigDI, &EMatR, &EMatI );
955  }
956  }
957  }
958 
959  //================================================ Report completion
960  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, "Rotation function distance computation complete." );
961 
962  //================================================ Done
963  return ( ret );
964 
965 }
ProSHADE_internal_distances::computeWeightsForEMatricesForLM
proshade_double computeWeightsForEMatricesForLM(ProSHADE_internal_data::ProSHADE_data *obj1, ProSHADE_internal_data::ProSHADE_data *obj2, proshade_unsign bandIter, proshade_unsign orderIter, proshade_double *obj1Vals, proshade_double *obj2Vals, proshade_unsign integOrder, proshade_double *abscissas, proshade_double *weights, proshade_single sphereDist)
This function computes the E matrix weight values for a given band and order and saves these into the...
Definition: ProSHADE_distances.cpp:433
ProSHADE_internal_distances::normaliseEMatrices
void normaliseEMatrices(ProSHADE_internal_data::ProSHADE_data *obj1, ProSHADE_internal_data::ProSHADE_data *obj2, ProSHADE_settings *settings)
This function normalises the E matrices.
Definition: ProSHADE_distances.cpp:572
ProSHADE_settings::integOrder
proshade_unsign integOrder
The order required for full Gauss-Legendre integration between the spheres.
Definition: ProSHADE_settings.hpp:68
ProSHADE_internal_data::ProSHADE_data::getEMatrixByBand
proshade_complex ** getEMatrixByBand(proshade_unsign band)
This function allows access to E matrix for a particular band.
Definition: ProSHADE_data.cpp:3827
ProSHADE_internal_maths::gaussLegendreIntegration
void gaussLegendreIntegration(proshade_complex *vals, proshade_unsign valsSize, proshade_unsign order, proshade_double *abscissas, proshade_double *weights, proshade_double integralOverRange, proshade_double maxSphereDists, proshade_double *retReal, proshade_double *retImag)
Function to compute the complete complex Gauss-Legendre integration over spherical harmonic values in...
Definition: ProSHADE_maths.cpp:711
ProSHADE_settings::computeTraceSigmaDesc
bool computeTraceSigmaDesc
If true, the trace sigma descriptor will be computed, otherwise all its computations will be omitted.
Definition: ProSHADE_settings.hpp:110
ProSHADE_distances.hpp
This is the header file containing declarations of functions required for computation of shape distan...
ProSHADE_internal_maths::gaussLegendreIntegrationReal
proshade_double gaussLegendreIntegrationReal(proshade_double *vals, proshade_unsign valsSize, proshade_unsign order, proshade_double *abscissas, proshade_double *weights, proshade_double integralOverRange, proshade_double maxSphereDists)
Function to compute real part of the Gauss-Legendre integration over spherical harmonic values in dif...
Definition: ProSHADE_maths.cpp:621
ProSHADE_settings::computeRotationFuncDesc
bool computeRotationFuncDesc
If true, the rotation function descriptor will be computed, otherwise all its computations will be om...
Definition: ProSHADE_settings.hpp:111
ProSHADE_settings::taylorSeriesCap
proshade_unsign taylorSeriesCap
The max limit on the Taylor series expansion done for the abscissas of the Gauss-Legendre integration...
Definition: ProSHADE_settings.hpp:69
ProSHADE_internal_data::ProSHADE_data::getComparisonBand
proshade_unsign getComparisonBand(void)
This function allows access to the maximum band for the comparison.
Definition: ProSHADE_data.cpp:3879
ProSHADE_internal_distances::allocateTrSigmaWorkspace
void allocateTrSigmaWorkspace(proshade_unsign minSpheres, proshade_unsign intOrder, proshade_double *&obj1Vals, proshade_double *&obj2Vals, proshade_double *&GLabscissas, proshade_double *&glWeights, proshade_complex *&radiiVals)
This helper function is responsible for allocating the workspace memory required for trace sigma desc...
Definition: ProSHADE_distances.cpp:312
ProSHADE_internal_data::ProSHADE_data::setIntegrationWeight
void setIntegrationWeight(proshade_double intW)
This function allows setting the integration weight for the object.
Definition: ProSHADE_data.cpp:4097
ProSHADE_internal_data::ProSHADE_data::normaliseEMatrixValue
void normaliseEMatrixValue(proshade_unsign band, proshade_unsign order1, proshade_unsign order2, proshade_double normF)
This function allows normalising the E matrix value.
Definition: ProSHADE_data.cpp:4146
ProSHADE_exception
This class is the representation of ProSHADE exception.
Definition: ProSHADE_exceptions.hpp:37
ProSHADE_internal_data::ProSHADE_data::allocateRRPMemory
void allocateRRPMemory()
This function allocates the required memory for the RRP matrices.
Definition: ProSHADE_distances.cpp:31
ProSHADE_internal_data::ProSHADE_data::getEMatrixValue
void getEMatrixValue(proshade_unsign band, proshade_unsign order1, proshade_unsign order2, proshade_double *valueReal, proshade_double *valueImag)
This function allows access to E matrix by knowing the band, order1 and order2 indices.
Definition: ProSHADE_data.cpp:3842
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::noSpheres
proshade_unsign noSpheres
The number of spheres with map projected onto them.
Definition: ProSHADE_data.hpp:119
ProSHADE_internal_distances::prepareInvSOFTPlan
void prepareInvSOFTPlan(fftw_plan *inverseSO3, int band, fftw_complex *work1, proshade_complex *invCoeffs)
This function prepares the FFTW plan for the inverse SO(3) transform.
Definition: ProSHADE_distances.cpp:785
ProSHADE_internal_distances::releaseInvSOFTMemory
void releaseInvSOFTMemory(proshade_complex *&work1, proshade_complex *&work2, proshade_double *&work3)
This function releases the memory used for computation of the inverse SOFT transform.
Definition: ProSHADE_distances.cpp:832
ProSHADE_internal_data::ProSHADE_data::getMaxSpheres
proshade_unsign getMaxSpheres(void)
This function returns the number of spheres which contain the whole object.
Definition: ProSHADE_data.cpp:3609
ProSHADE_internal_data::ProSHADE_data::setSO3CoeffValue
void setSO3CoeffValue(proshade_unsign position, proshade_complex val)
This function allows setting the SOFT coefficient values using array position and value.
Definition: ProSHADE_data.cpp:4162
ProSHADE_internal_distances::isBandWithinShell
bool isBandWithinShell(proshade_unsign bandInQuestion, proshade_unsign shellInQuestion, ProSHADE_internal_spheres::ProSHADE_sphere **spheres)
This function checks if a band is available for a given shell.
Definition: ProSHADE_distances.cpp:139
ProSHADE_internal_data::ProSHADE_data::computeRRPMatrices
void computeRRPMatrices(ProSHADE_settings *settings)
This function pre-computes the RRP matrices for a data object.
Definition: ProSHADE_distances.cpp:59
ProSHADE_internal_data::ProSHADE_data::getAnySphereRadius
proshade_double getAnySphereRadius(proshade_unsign shell)
This function allows access to the radius of any particular sphere.
Definition: ProSHADE_data.cpp:3780
ProSHADE_settings::progressiveSphereMapping
bool progressiveSphereMapping
If true, each shell will have its own angular resolution dependent on the actual number of map points...
Definition: ProSHADE_settings.hpp:102
ProSHADE_internal_data::ProSHADE_data::getSO3Coeffs
proshade_complex * getSO3Coeffs(void)
This function allows access to the SO(3) coefficients array.
Definition: ProSHADE_data.cpp:3868
ProSHADE_settings::maxSphereDists
proshade_single maxSphereDists
The distance between spheres in spherical mapping for the largest sphere.
Definition: ProSHADE_settings.hpp:65
ProSHADE_internal_data::ProSHADE_data::maxShellBand
proshade_unsign maxShellBand
The maximum band for any shell of the object.
Definition: ProSHADE_data.hpp:123
ProSHADE_internal_data::ProSHADE_data::getIntegrationWeight
proshade_double getIntegrationWeight(void)
This function allows access to the integration weight for the object.
Definition: ProSHADE_data.cpp:3791
ProSHADE_internal_data::ProSHADE_data::setIntegrationWeightCumul
void setIntegrationWeightCumul(proshade_double intW)
This function allows setting the cumulative integration weight for the object.
Definition: ProSHADE_data.cpp:4111
ProSHADE_internal_spheres::ProSHADE_sphere
This class contains all inputed and derived data for a single sphere.
Definition: ProSHADE_spheres.hpp:49
ProSHADE_internal_distances::computeRotationunctionDescriptor
proshade_double computeRotationunctionDescriptor(ProSHADE_internal_data::ProSHADE_data *obj1, ProSHADE_internal_data::ProSHADE_data *obj2, ProSHADE_settings *settings)
This function computes the rotation function descriptor value between two objects.
Definition: ProSHADE_distances.cpp:906
ProSHADE_settings::verbose
proshade_signed verbose
Should the software report on the progress, or just be quiet? Value between -1 (nothing) and 4 (loud)
Definition: ProSHADE_settings.hpp:140
ProSHADE_internal_maths::complexMatrixSVDSigmasOnly
void complexMatrixSVDSigmasOnly(proshade_complex **mat, int dim, double *&singularValues)
Function to compute the complete complex matrix SVD and return only the sigmas.
Definition: ProSHADE_maths.cpp:804
ProSHADE_internal_data::ProSHADE_data::getWignerMatrixValue
void getWignerMatrixValue(proshade_unsign band, proshade_unsign order1, proshade_unsign order2, proshade_double *valueReal, proshade_double *valueImag)
This function allows access to the Wigner D matrix by knowing the band, order1 and order2 indices.
Definition: ProSHADE_data.cpp:3894
ProSHADE_internal_data::ProSHADE_data::getShellBandwidth
proshade_unsign getShellBandwidth(proshade_unsign shell)
This function allows access to the bandwidth of a particular shell.
Definition: ProSHADE_data.cpp:3803
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_distances::computeRRPPearsonCoefficients
void computeRRPPearsonCoefficients(ProSHADE_internal_data::ProSHADE_data *obj1, ProSHADE_internal_data::ProSHADE_data *obj2, ProSHADE_settings *settings, proshade_unsign minCommonBands, proshade_unsign minCommonShells, std::vector< proshade_double > *bandDists)
This function gets the Pearson's coefficients or all bands between two objects.
Definition: ProSHADE_distances.cpp:211
ProSHADE_internal_maths::pearsonCorrCoeff
proshade_double pearsonCorrCoeff(proshade_double *valSet1, proshade_double *valSet2, proshade_unsign length)
Function for computing the Pearson's correlation coefficient.
Definition: ProSHADE_maths.cpp:246
ProSHADE_internal_wigner::computeWignerMatricesForRotation
void computeWignerMatricesForRotation(ProSHADE_settings *settings, ProSHADE_internal_data::ProSHADE_data *obj, proshade_double eulerAlpha, proshade_double eulerBeta, proshade_double eulerGamma)
This function computes the Wigner D matrices for a particular set of Euler angles.
Definition: ProSHADE_wignerMatrices.cpp:258
ProSHADE_internal_data::ProSHADE_data::rrpMatrices
proshade_double *** rrpMatrices
The energy levels descriptor shell correlation tables.
Definition: ProSHADE_data.hpp:126
ProSHADE_internal_peakSearch::getBestPeakEulerAngsNaive
void getBestPeakEulerAngsNaive(proshade_complex *map, proshade_unsign dim, proshade_double *eulA, proshade_double *eulB, proshade_double *eulG, ProSHADE_settings *settings)
This function finds the highest peaks optimised Euler angles using the "naive" approach.
Definition: ProSHADE_peakSearch.cpp:352
ProSHADE_settings::task
ProSHADE_Task task
This custom type variable determines which task to perfom (i.e. symmetry detection,...
Definition: ProSHADE_settings.hpp:40
ProSHADE_internal_data::ProSHADE_data::getRealSphHarmValue
proshade_double * getRealSphHarmValue(proshade_unsign band, proshade_unsign order, proshade_unsign shell)
This function allows access to the private internal real spherical harmonics values.
Definition: ProSHADE_data.cpp:3754
ProSHADE_internal_maths::complexMultiplicationConjugRealOnly
proshade_double complexMultiplicationConjugRealOnly(proshade_double *r1, proshade_double *i1, proshade_double *r2, proshade_double *i2)
Function to conjuggate multiply two complex numbers and return the real part only.
Definition: ProSHADE_maths.cpp:103
ProSHADE_internal_data::ProSHADE_data::shellBandExists
bool shellBandExists(proshade_unsign shell, proshade_unsign bandVal)
This function checks if particular shell has a particular band.
Definition: ProSHADE_data.cpp:3656
ProSHADE_internal_data::ProSHADE_data::getMaxBand
proshade_unsign getMaxBand(void)
This function returns the maximum band value for the object.
Definition: ProSHADE_data.cpp:3630
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::getImagSphHarmValue
proshade_double * getImagSphHarmValue(proshade_unsign band, proshade_unsign order, proshade_unsign shell)
This function allows access to the private internal imaginary spherical harmonics values.
Definition: ProSHADE_data.cpp:3767
ProSHADE_internal_distances::computeEMatrices
void computeEMatrices(ProSHADE_internal_data::ProSHADE_data *obj1, ProSHADE_internal_data::ProSHADE_data *obj2, ProSHADE_settings *settings)
This function computes the complete E matrices and their weights between any two objects.
Definition: ProSHADE_distances.cpp:510
ProSHADE_internal_data::ProSHADE_data::getSpherePosValue
proshade_single getSpherePosValue(proshade_unsign shell)
This function allows access to sphere positions.
Definition: ProSHADE_data.cpp:3815
ProSHADE_internal_data::ProSHADE_data::getRRPValue
proshade_double getRRPValue(proshade_unsign band, proshade_unsign sh1, proshade_unsign sh2)
This function allows access to the priva internal RRP matrices.
Definition: ProSHADE_data.cpp:3640
ProSHADE_internal_data::ProSHADE_data::allocateEMatrices
void allocateEMatrices(proshade_unsign band)
This function allocates the required memory for the E matrices.
Definition: ProSHADE_distances.cpp:275
ProSHADE_internal_distances::computeInverseSOFTTransform
void computeInverseSOFTTransform(ProSHADE_internal_data::ProSHADE_data *obj1, ProSHADE_internal_data::ProSHADE_data *obj2, ProSHADE_settings *settings)
This function computes the inverse SO(3) transform.
Definition: ProSHADE_distances.cpp:854
ProSHADE_internal_distances::computeEMatricesForLM
void computeEMatricesForLM(ProSHADE_internal_data::ProSHADE_data *obj1, ProSHADE_internal_data::ProSHADE_data *obj2, proshade_unsign bandIter, proshade_unsign orderIter, proshade_complex *radiiVals, proshade_unsign integOrder, proshade_double *abscissas, proshade_double *weights, proshade_double integRange, proshade_double sphereDist)
This function computes the E matrix un-weighted values for a given band and order and saves these int...
Definition: ProSHADE_distances.cpp:370
ProSHADE_internal_maths::getLegendreAbscAndWeights
void getLegendreAbscAndWeights(proshade_unsign order, proshade_double *abscissas, proshade_double *weights, proshade_unsign taylorSeriesCap)
Function to prepare abscissas and weights for Gauss-Legendre integration.
Definition: ProSHADE_maths.cpp:289
ProSHADE_internal_distances::computeSphericalHarmonicsMagnitude
void computeSphericalHarmonicsMagnitude(ProSHADE_internal_data::ProSHADE_data *obj, proshade_unsign band, proshade_unsign order, proshade_unsign radius, proshade_double *result)
This function computes the magnitude of a particular spherical harmonics position for a given object,...
Definition: ProSHADE_distances.cpp:341
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_misc::checkMemoryAllocation
void checkMemoryAllocation(chVar checkVar, std::string fileP, unsigned int lineP, std::string funcP, std::string infoP="This error may occurs when ProSHADE requests memory to be\n : allocated to it and this operation fails. This could\n : happen when not enough memory is available, either due to\n : other processes using a lot of memory, or when the machine\n : does not have sufficient memory available. Re-run to see\n : if this problem persists.")
Checks if memory was allocated properly.
Definition: ProSHADE_misc.hpp:67
ProSHADE_internal_distances::releaseTrSigmaWorkspace
void releaseTrSigmaWorkspace(proshade_double *&obj1Vals, proshade_double *&obj2Vals, proshade_double *&GLabscissas, proshade_double *&glWeights, proshade_complex *&radiiVals)
This helper function is responsible for deleting the workspace memory required for trace sigma descri...
Definition: ProSHADE_distances.cpp:478
ProSHADE_internal_distances::allocateInvSOFTWorkspaces
void allocateInvSOFTWorkspaces(proshade_complex *&work1, proshade_complex *&work2, proshade_double *&work3, proshade_unsign band)
This function allocates the workspaces required to compute the inverse SOFT transform.
Definition: ProSHADE_distances.cpp:761
ProSHADE_internal_data::ProSHADE_data::getInvSO3Coeffs
proshade_complex * getInvSO3Coeffs(void)
This function allows access to the inverse SO(3) coefficients array.
Definition: ProSHADE_data.cpp:3857
ProSHADE_internal_maths::complexMultiplicationConjug
void complexMultiplicationConjug(proshade_double *r1, proshade_double *i1, proshade_double *r2, proshade_double *i2, proshade_double *retReal, proshade_double *retImag)
Function to multiply two complex numbers by using the second number's conjugate.
Definition: ProSHADE_maths.cpp:62
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_distances::generateSO3CoeffsFromEMatrices
void generateSO3CoeffsFromEMatrices(ProSHADE_internal_data::ProSHADE_data *obj1, ProSHADE_internal_data::ProSHADE_data *obj2, ProSHADE_settings *settings)
This function converts the E matrices to SO(3) coefficients.
Definition: ProSHADE_distances.cpp:701
ProSHADE_internal_data::ProSHADE_data::setEMatrixValue
void setEMatrixValue(proshade_unsign band, proshade_unsign order1, proshade_unsign order2, proshade_complex val)
This function allows setting the E matrix value.
Definition: ProSHADE_data.cpp:4128
ProSHADE_internal_data::ProSHADE_data::allocateSO3CoeffsSpace
void allocateSO3CoeffsSpace(proshade_unsign band)
This function allocates the memory for the SO(3) coefficients and the inverse for that calling object...
Definition: ProSHADE_distances.cpp:675
ProSHADE_internal_messages::printProgressMessage
void printProgressMessage(proshade_signed verbose, proshade_signed messageLevel, std::string message)
General stdout message printing.
Definition: ProSHADE_messages.cpp:70
ProSHADE_settings::enLevMatrixPowerWeight
proshade_double enLevMatrixPowerWeight
If RRP matrices shell position is to be weighted by putting the position as an exponent,...
Definition: ProSHADE_settings.hpp:109
ProSHADE_internal_maths::complexMultiplicationRealOnly
proshade_double complexMultiplicationRealOnly(proshade_double *r1, proshade_double *i1, proshade_double *r2, proshade_double *i2)
Function to multiply two complex numbers and return the real part only.
Definition: ProSHADE_maths.cpp:83
ProSHADE_settings::computeEnergyLevelsDesc
bool computeEnergyLevelsDesc
If true, the energy levels descriptor will be computed, otherwise all its computations will be omitte...
Definition: ProSHADE_settings.hpp:108