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

This namespace contains the functions used for computing distances between two structures as represented by the ProSHADE_data class. More...

Functions

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. More...
 
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. More...
 
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. More...
 
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 descriptor computation. More...
 
void computeSphericalHarmonicsMagnitude (ProSHADE_internal_data::ProSHADE_data *obj, int band, int order, proshade_unsign radius, proshade_double *result)
 This function computes the magnitude of a particular spherical harmonics position for a given object, weighting it by the radius^2 (for integration). More...
 
void computeEMatricesForLM (ProSHADE_internal_data::ProSHADE_data *obj1, ProSHADE_internal_data::ProSHADE_data *obj2, int bandIter, int orderIter, proshade_complex *radiiVals, int 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 into the obj2 parameter. More...
 
proshade_double computeWeightsForEMatricesForLM (ProSHADE_internal_data::ProSHADE_data *obj1, ProSHADE_internal_data::ProSHADE_data *obj2, int bandIter, int orderIter, proshade_double *obj1Vals, proshade_double *obj2Vals, int 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 appropriate objects. More...
 
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 descriptor computation. More...
 
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. More...
 
void normaliseEMatrices (ProSHADE_internal_data::ProSHADE_data *obj1, ProSHADE_internal_data::ProSHADE_data *obj2, ProSHADE_settings *settings)
 This function normalises the E matrices. More...
 
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. More...
 
void generateSO3CoeffsFromEMatrices (ProSHADE_internal_data::ProSHADE_data *obj2, ProSHADE_settings *settings)
 This function converts the E matrices to SO(3) coefficients. More...
 
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. More...
 
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. More...
 
void releaseInvSOFTMemory (proshade_complex *&work1, proshade_complex *&work2, proshade_double *&work3)
 This function releases the memory used for computation of the inverse SOFT transform. More...
 
void computeInverseSOFTTransform (ProSHADE_internal_data::ProSHADE_data *obj2, ProSHADE_settings *settings)
 This function computes the inverse SO(3) transform. More...
 
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. More...
 

Detailed Description

This namespace contains the functions used for computing distances between two structures as represented by the ProSHADE_data class.

The ProSHADE_internal_distances namespace contains all the functionality required to compute distances between two or more ProSHADE_data objects and also all of the pre-computation functions.

Function Documentation

◆ allocateInvSOFTWorkspaces()

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

Parameters
[in]work1The first workspace pointer to be allocated.
[in]work2The second workspace pointer to be allocated.
[in]work3The third workspace pointer to be allocated.
[in]bandThe bandwidth of the computations (this determines how much workspace will be required).

Definition at line 792 of file ProSHADE_distances.cpp.

793 {
794  //================================================ Allocate memory
795  work1 = new proshade_complex[8 * static_cast<proshade_unsign> ( pow( static_cast<double> ( band ), 3.0 ) )];
796  work2 = new proshade_complex[14 * static_cast<proshade_unsign> ( pow( static_cast<double> ( band ), 2.0 ) ) + (48 * band)];
797  work3 = new proshade_double [2 * static_cast<proshade_unsign> ( pow( static_cast<double> ( band ), 2.0 ) ) + (24 * band)];
798 
799  //================================================ Check the memory allocation
800  ProSHADE_internal_misc::checkMemoryAllocation ( work1, __FILE__, __LINE__, __func__ );
801  ProSHADE_internal_misc::checkMemoryAllocation ( work2, __FILE__, __LINE__, __func__ );
802  ProSHADE_internal_misc::checkMemoryAllocation ( work3, __FILE__, __LINE__, __func__ );
803 
804  //================================================ Done
805  return ;
806 
807 }

◆ allocateTrSigmaWorkspace()

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 
)

This helper function is responsible for allocating the workspace memory required for trace sigma descriptor computation.

Parameters
[in]minSpheresThe minima of the number of spheres available in the compared objects.
[in]intOrderThe integration order for the computation.
[in]obj1ValsArray to hold the shell values for the first object integgration.
[in]obj2ValsArray to hold the shell values for the second object integgration.
[in]GLabscissasAn array to hold the pre-computed anscissas for the Gauss-Legendre integration.
[in]glWeightsAn array to hold the pre-computed weights for the Gauss-Legendre integration.
[in]radiiValsA complex array to hold the results of combining spherical harmonics coefficients of the two objects for each shell.

Definition at line 329 of file ProSHADE_distances.cpp.

330 {
331  //================================================ Allocate the memory
332  obj1Vals = new proshade_double [minSpheres];
333  obj2Vals = new proshade_double [minSpheres];
334  radiiVals = new proshade_complex[minSpheres];
335  GLabscissas = new proshade_double [intOrder];
336  GLweights = new proshade_double [intOrder];
337 
338  //================================================ Check the memory allocation
339  ProSHADE_internal_misc::checkMemoryAllocation ( obj1Vals, __FILE__, __LINE__, __func__ );
340  ProSHADE_internal_misc::checkMemoryAllocation ( obj2Vals, __FILE__, __LINE__, __func__ );
341  ProSHADE_internal_misc::checkMemoryAllocation ( radiiVals, __FILE__, __LINE__, __func__ );
342  ProSHADE_internal_misc::checkMemoryAllocation ( GLabscissas, __FILE__, __LINE__, __func__ );
343  ProSHADE_internal_misc::checkMemoryAllocation ( GLweights, __FILE__, __LINE__, __func__ );
344 
345  //================================================ Done
346  return ;
347 
348 }

◆ computeEMatrices()

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

This function allocates the space required for storing the E matrices, allocates all the workspace requierd for the computation and proceeds to compute the values for all band (l), order1(m) and order2(m') E matrix values. It then proceeds to release all non required memory and terminates, leaving all its results in the second ProSHADE data object supplied. This function does NOT apply the weights to the matrices, it needs to be done subsequently!

Parameters
[in]obj1The first ProSHADE_data object for which the computation is done.
[in]obj2The second ProSHADE_data object for which the computation is done.
[in]settingsA pointer to settings class containing all the information required for the task.

Definition at line 540 of file ProSHADE_distances.cpp.

541 {
542  //================================================ Report progress
543  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, "Starting computation of E matrices.", settings->messageShift );
544 
545  //================================================ Allocatre memory for E matrices in the second object (first may be compared to more structures and therefore its data would be written over)
546  obj2->allocateEMatrices ( std::min ( obj1->getMaxBand(), obj2->getMaxBand() ), settings->resolutionOversampling );
547 
548  //================================================ Initialise local variables
549  proshade_double *obj1Vals, *obj2Vals, *GLAbscissas, *GLWeights;
550  proshade_complex* radiiVals;
551  proshade_double integRange;
552 
553  //================================================ Allocate workspace memory
554  allocateTrSigmaWorkspace ( std::min( obj1->getMaxSpheres(), obj2->getMaxSpheres() ), settings->integOrder, obj1Vals, obj2Vals, GLAbscissas, GLWeights, radiiVals);
555 
556  //================================================ For each band (l), compute the E matrix integrals
557  for ( int bandIter = 0; bandIter < static_cast< int > ( std::min ( obj1->getMaxBand(), obj2->getMaxBand() ) ); bandIter++ )
558  {
559  //============================================ Initialise abscissas and weights for integration
560  int localIntegOrder = std::min ( static_cast< int > ( std::max ( bandIter, 5 ) ), static_cast< int > ( settings->integOrder ) );
561  if ( settings->noIntegrationSpeedup ) { localIntegOrder = static_cast< int > ( settings->integOrder ); }
562  ProSHADE_internal_maths::getLegendreAbscAndWeights ( static_cast< proshade_unsign > ( localIntegOrder ), GLAbscissas, GLWeights, settings->integApproxSteps );
563 
564  //============================================ For each order (m)
565  for ( int orderIter = 0; orderIter < ( ( bandIter * 2 ) + 1 ); orderIter++ )
566  {
567  //======================================== Get weights for the required band(l) and order (m)
568  integRange = computeWeightsForEMatricesForLM ( obj1, obj2, bandIter, orderIter, obj1Vals, obj2Vals, localIntegOrder, GLAbscissas, GLWeights, settings->maxSphereDists );
569 
570  //======================================== Compute E matrices value for given band (l) and order(m)
571  computeEMatricesForLM ( obj1, obj2, bandIter, orderIter, radiiVals, localIntegOrder, GLAbscissas, GLWeights, integRange, static_cast< proshade_double > ( settings->maxSphereDists ) );
572  }
573 
574  //============================================ Report progress
575  if ( settings->verbose > 3 )
576  {
577  std::stringstream hlpSS;
578  hlpSS << "E matrices computed for band " << bandIter;
579  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 4, hlpSS.str(), settings->messageShift );
580  }
581  }
582 
583  //================================================ Release the workspace memory
584  releaseTrSigmaWorkspace ( obj1Vals, obj2Vals, GLAbscissas, GLWeights, radiiVals );
585 
586  //================================================ Report progress
587  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 3, "E matrices computed.", settings->messageShift );
588 
589  //================================================ Done
590  return ;
591 
592 }

◆ computeEMatricesForLM()

void ProSHADE_internal_distances::computeEMatricesForLM ( ProSHADE_internal_data::ProSHADE_data obj1,
ProSHADE_internal_data::ProSHADE_data obj2,
int  bandIter,
int  orderIter,
proshade_complex *  radiiVals,
int  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 into the obj2 parameter.

Parameters
[in]obj1The ProSHADE_data object for which the comparison is done in regards to.
[in]obj2The ProSHADE_data object for which the comparison is done in regards from - the E matrices will be saved into this object.
[in]bandThe bandwidth of the SH value for which this should be done.
[in]orderThe order of the SH value for which this should be done.
[in]radiiValsAlready allocated array of proshade_complex to which integrated values will be saved. It must have size equal to minimum of spheres in the two compared objects.
[in]integOrderThe Gauss-Legendre integration order to be used.
[in]abscissasThe pre-computed abscissas for the Gauss-Legendre integration.
[in]weightsThe pre-computed weights for the Gauss-Legendre integration.
[in]integRangeThe range in angstroms between the smalleds and largest shell which are integrated over (might not be 0 to max for progressive shell sampling).
[in]sphereDistThe distance between any two spheres.

Definition at line 391 of file ProSHADE_distances.cpp.

392 {
393  //================================================ Initialise local variables
394  proshade_unsign objCombValsIter = 0;
395  proshade_double hlpReal, hlpImag, rSquared;
396  proshade_complex arrVal;
397  int o1ArrPos, o2ArrPos, locBand;
398  proshade_unsign integOrderU = static_cast< proshade_unsign > ( integOrder );
399 
400  //================================================ For each combination of m and m' for E matrices
401  for ( int order2Iter = 0; order2Iter < ( ( bandIter * 2 ) + 1 ); order2Iter++ )
402  {
403  //============================================ Reset loop
404  objCombValsIter = 0;
405 
406  //============================================ Find the c*conj(c) values for different radii
407  for ( proshade_unsign radiusIter = 0; radiusIter < std::min( obj1->getMaxSpheres(), obj2->getMaxSpheres() ); radiusIter++ )
408  {
409  //======================================== Get only values where the shell has the band
410  if ( std::min ( obj1->getShellBandwidth ( radiusIter ), obj2->getShellBandwidth ( radiusIter ) ) <= static_cast< proshade_unsign > ( bandIter ) ) { continue; }
411 
412  //======================================== Pre-compute values
413  rSquared = pow ( ( static_cast<proshade_double> ( obj1->getAnySphereRadius( radiusIter ) ) ), 2.0 );
414 
415  //======================================== Multiply coeffs
416  locBand = static_cast< int > ( obj1->spheres[radiusIter]->getLocalBandwidth() );
417  o1ArrPos = seanindex ( orderIter - bandIter, bandIter, locBand );
418  o2ArrPos = seanindex ( order2Iter - bandIter, bandIter, locBand );
419 
421  &obj1->sphericalHarmonics[radiusIter][o1ArrPos][1],
422  &obj2->sphericalHarmonics[radiusIter][o2ArrPos][0],
423  &obj2->sphericalHarmonics[radiusIter][o2ArrPos][1],
424  &hlpReal, &hlpImag );
425 
426  //======================================== Apply r^2 integral weight
427  radiiVals[objCombValsIter][0] = hlpReal * rSquared;
428  radiiVals[objCombValsIter][1] = hlpImag * rSquared;
429 
430  objCombValsIter += 1;
431  }
432 
433  //============================================ Integrate over all radii using n-point Gauss-Legendre integration
434  ProSHADE_internal_maths::gaussLegendreIntegration ( radiiVals, objCombValsIter, integOrderU, abscissas, weights, integRange, sphereDist, &hlpReal, &hlpImag );
435 
436  //============================================ Save the result into E matrices
437  arrVal[0] = hlpReal;
438  arrVal[1] = hlpImag;
439  obj2->setEMatrixValue ( bandIter, orderIter, order2Iter, arrVal );
440  }
441 
442  //================================================ Done
443  return ;
444 
445 }

◆ computeEnergyLevelsDescriptor()

proshade_double ProSHADE_internal_distances::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.

This function is where the enery levels descriptor computation is controlled and done from. It starts by making sure that both input data objects have the RRP matrices computed and then it proceeds to compute the Pearson's coefficients for each band, finally averaging the band values and returning the descriptos.

Parameters
[in]obj1The first ProSHADE_data object against which comparison is done.
[in]obj2The second ProSHADE_data object which is compared to the first.
[in]settingsA pointer to settings class containing all the information required for the task.

Definition at line 161 of file ProSHADE_distances.cpp.

162 {
163  //================================================ Report starting the task
164  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 1, "Starting energy levels distance computation.", settings->messageShift );
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.", settings->messageShift );
193 
194  //================================================ Done
195  return ( ret );
196 
197 }

◆ computeInverseSOFTTransform()

void ProSHADE_internal_distances::computeInverseSOFTTransform ( ProSHADE_internal_data::ProSHADE_data obj2,
ProSHADE_settings settings 
)

This function computes the inverse SO(3) transform.

This function firstly allocates all the required workspaces for the inverse SO(3) Fourier Transform, then it prepares the FFTW plans for performing the FFTW inverse Fourier transform in the SO(3) space using FFTW and finally it subjects the SO(3) coeffficients available at this point to the computation. The results are saved into the second object, memory is released and function terminates.

Parameters
[in]obj1The first ProSHADE_data object against which comparison is done.
[in]obj2The second ProSHADE_data object which is compared to the first.
[in]settingsA pointer to settings class containing all the information required for the task.

Definition at line 885 of file ProSHADE_distances.cpp.

886 {
887  //================================================ Report progress
888  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, "Computing inverse SO(3) Fourier transform.", settings->messageShift );
889 
890  //================================================ Initialise local variables
891  proshade_complex *workspace1, *workspace2;
892  proshade_double *workspace3;
893  fftw_plan inverseSO3;
894 
895  //================================================ Allocate memory for the workspaces
896  allocateInvSOFTWorkspaces ( workspace1, workspace2, workspace3, obj2->getEMatDim ( ) );
897 
898  //================================================ Prepare the FFTW plan
899  prepareInvSOFTPlan ( &inverseSO3, static_cast< int > ( obj2->getEMatDim ( ) ), workspace1, obj2->getInvSO3Coeffs ( ) );
900 
901  //================================================ Compute the transform
902  Inverse_SO3_Naive_fftw ( static_cast< int > ( obj2->getEMatDim ( ) ),
903  obj2->getSO3Coeffs ( ),
904  obj2->getInvSO3Coeffs ( ),
905  workspace1,
906  workspace2,
907  workspace3,
908  &inverseSO3,
909  0 );
910 
911  //================================================ Release memory
912  releaseInvSOFTMemory ( workspace1, workspace2, workspace3 );
913  fftw_destroy_plan ( inverseSO3 );
914 
915  //================================================ Report progress
916  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 3, "Inverse SO(3) Fourier transform computed.", settings->messageShift );
917 
918 
919 
920  //================================================ Done
921  return ;
922 
923 }

◆ computeRotationFunctionDescriptor()

proshade_double ProSHADE_internal_distances::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.

This function starts by sanity checks and computation and normalisation of the E matrices used for the trace sigma descriptor (these only need to be computed if trace sigma descriptor is NOT computed, otherwise the stored values are used to save time). It then converts the E matrix values to SO(3) transform coefficients and invers these using the SOFT library. From the resulting rotation function map, it selects the highest peak and applies its rotation to the E matrix values. The resulting values are simply summed, as this sum can be proven to be the argmin of the distance between the two objects in their spherical harmonics decomposition space.

Parameters
[in]obj1The first ProSHADE_data object against which comparison is done.
[in]obj2The second ProSHADE_data object which is compared to the first.
[in]settingsA pointer to settings class containing all the information required for the task.
[out]retThe final normalised value of the rotation function descriptor for the two objects.

Definition at line 939 of file ProSHADE_distances.cpp.

940 {
941  //================================================ Report starting the task
942  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 1, "Starting rotation function distance computation.", settings->messageShift );
943 
944  //================================================ Initialise return variable
945  proshade_double ret = 0.0;
946  proshade_double eulA, eulB, eulG, EMatR, EMatI, WigDR, WigDI;
947 
948  //================================================ Sanity check
949  if ( !settings->computeRotationFuncDesc )
950  {
951  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!" );
952  }
953 
954  //================================================ Compute weighted E matrices if not already present
955  if ( !settings->computeTraceSigmaDesc )
956  {
957  computeEMatrices ( obj1, obj2, settings );
958  normaliseEMatrices ( obj1, obj2, settings );
959  }
960 
961  //================================================ Generate SO(3) coefficients
962  generateSO3CoeffsFromEMatrices ( obj2, settings );
963 
964  //================================================ Compute the inverse SO(3) Fourier Transform (SOFT) on the newly computed coefficients
965  computeInverseSOFTTransform ( obj2, settings );
966 
967  //================================================ Get inverse SO(3) map top peak Euler angle values
969  obj2->getEMatDim ( ) * 2,
970  &eulA, &eulB, &eulG, settings );
971 
972  //================================================ Compute the Wigner D matrices for the Euler angles
973  ProSHADE_internal_wigner::computeWignerMatricesForRotation ( settings, obj2, eulA, eulB, eulG );
974 
975  //================================================ Compute the distance
976  for ( proshade_unsign bandIter = 0; bandIter < obj2->getEMatDim(); bandIter++ )
977  {
978  //============================================ For each order1
979  for ( proshade_unsign order1 = 0; order1 < ( ( bandIter * 2 ) + 1 ); order1++ )
980  {
981  //======================================== For each order2
982  for ( proshade_unsign order2 = 0; order2 < ( ( bandIter * 2 ) + 1 ); order2++ )
983  {
984  //==================================== Multiply D_{l} * E_{l} and get sum over l of traces (i.e. just sum all together)
985  obj2->getEMatrixValue ( bandIter, order1, order2, &EMatR, &EMatI );
986  obj2->getWignerMatrixValue ( bandIter, order2, order1, &WigDR, &WigDI );
987  ret += ProSHADE_internal_maths::complexMultiplicationRealOnly ( &WigDR, &WigDI, &EMatR, &EMatI );
988  }
989  }
990  }
991 
992  //================================================ Report completion
993  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, "Rotation function distance computation complete.", settings->messageShift );
994 
995  //================================================ Done
996  return ( ret );
997 
998 }

◆ computeRRPPearsonCoefficients()

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 
)

This function gets the Pearson's coefficients or all bands between two objects.

This function takes two data objects with their RRP matrices computed and proceeds to compute the Pearson's correlation coefficient for each band, saving it into the supplied vector.

Parameters
[in]obj1The first ProSHADE_data object against which comparison is done.
[in]obj2The second ProSHADE_data object which is compared to the first.
[in]settingsA pointer to settings class containing all the information required for the task.
[in]minCommonBandsThe number of common bands between the two objects.
[in]minCommonShellsThe index of highest common shell in both objects.
[in]bandDistsEmpty vector of proshade_doubles to which the Pearson's Coefficients will be saved for each band.

Definition at line 211 of file ProSHADE_distances.cpp.

212 {
213  //================================================ Report completion
214  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, "Correlating RRP matrices.", settings->messageShift );
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.", settings->messageShift );
262 
263  //================================================ Done
264  return ;
265 }

◆ computeSphericalHarmonicsMagnitude()

void ProSHADE_internal_distances::computeSphericalHarmonicsMagnitude ( ProSHADE_internal_data::ProSHADE_data obj,
int  band,
int  order,
proshade_unsign  radius,
proshade_double *  result 
)

This function computes the magnitude of a particular spherical harmonics position for a given object, weighting it by the radius^2 (for integration).

Parameters
[in]objThe ProSHADE_data object for which the computation is to be done.
[in]bandThe bandwidth of the SH value for which this should be done.
[in]orderThe order of the SH value for which this should be done.
[in]radiusThe shell of the SH value for which this should be done.
[in]resultThe location where the result is to be saved.

Definition at line 358 of file ProSHADE_distances.cpp.

359 {
360  //================================================ Pre-compute values
361  int locBand = static_cast< int > ( obj->spheres[radius]->getLocalBandwidth() );
362  int objArrPos = seanindex ( order - band, band, locBand );
363 
364  //================================================ Find the magnitude
366  &obj->sphericalHarmonics[radius][objArrPos][1],
367  &obj->sphericalHarmonics[radius][objArrPos][0],
368  &obj->sphericalHarmonics[radius][objArrPos][1] );
369 
370  //================================================ Weight by radius^2 for the integration that will follow
371  *result *= pow ( static_cast<proshade_double> ( obj->getAnySphereRadius( radius ) ), 2.0 );
372 
373  //================================================ Done
374  return ;
375 
376 }

◆ computeTraceSigmaDescriptor()

proshade_double ProSHADE_internal_distances::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.

This function starts by checking if the trace sigma descriptor was requested and if so, proceeds to compute the E matrices. These are 3D matrices with each l,m,m' value being the combination of the c_{l,m} and c*_{l,m'} spherical harmonics coefficients. Once computed, the E matrices are normalised by the magnitudes of the objects spherical harmonics coefficients and the SVD is computed for each l (i.e. on each m x m' matrix). The sum of the trace of the sigmas of the SVD is then the trace sigma descriptor, whose value is returned.

Parameters
[in]obj1The first ProSHADE_data object against which comparison is done.
[in]obj2The second ProSHADE_data object which is compared to the first.
[in]settingsA pointer to settings class containing all the information required for the task.
[out]retThe final normalised value of the trace sigma descriptor for the two objects.

Definition at line 648 of file ProSHADE_distances.cpp.

649 {
650  //================================================ Report starting the task
651  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 1, "Starting trace sigma distance computation.", settings->messageShift );
652 
653  //================================================ Initialise return variable
654  proshade_double ret = 0.0;
655 
656  //================================================ Sanity check
657  if ( !settings->computeTraceSigmaDesc )
658  {
659  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!" );
660  }
661 
662  //================================================ Empty the cumulative weights back to 0.0 for each structure
663  obj1->setIntegrationWeight ( 0.0 );
664  obj1->setIntegrationWeight ( 0.0 );
665 
666  //================================================ Compute un-weighted E matrices and their weights
667  computeEMatrices ( obj1, obj2, settings );
668 
669  //================================================ Normalise E matrices by the magnitudes
670  normaliseEMatrices ( obj1, obj2, settings );
671 
672  //================================================ Allocate the required memory
673  double* singularValues = new double[( ( std::min ( obj1->getMaxBand(), obj2->getMaxBand() ) * 2 ) + 1 )];
674  ProSHADE_internal_misc::checkMemoryAllocation ( singularValues, __FILE__, __LINE__, __func__ );
675 
676  //================================================ Compute the distance
677  for ( proshade_unsign lIter = 0; lIter < std::min ( obj1->getMaxBand(), obj2->getMaxBand() ); lIter++ )
678  {
679  //============================================ Find the complex matrix SVD singular values
680  ProSHADE_internal_maths::complexMatrixSVDSigmasOnly ( obj2->getEMatrixByBand ( lIter ), static_cast<int> ( ( lIter * 2 ) + 1 ), singularValues );
681 
682  //============================================ Now sum the trace
683  for ( proshade_unsign iter = 0; iter < ( ( lIter * 2 ) + 1 ); iter++ )
684  {
685  ret += singularValues[iter];
686  }
687  }
688 
689  //================================================ Report completion
690  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 3, "E matrices decomposed to singular values.", settings->messageShift );
691 
692  //================================================ Release the memory
693  delete[] singularValues;
694 
695  //================================================ Report completion
696  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, "Trace sigma distance computation complete.", settings->messageShift );
697 
698  //================================================ Done
699  return ( ret );
700 
701 }

◆ computeWeightsForEMatricesForLM()

proshade_double ProSHADE_internal_distances::computeWeightsForEMatricesForLM ( ProSHADE_internal_data::ProSHADE_data obj1,
ProSHADE_internal_data::ProSHADE_data obj2,
int  bandIter,
int  orderIter,
proshade_double *  obj1Vals,
proshade_double *  obj2Vals,
int  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 appropriate objects.

Parameters
[in]obj1The ProSHADE_data object for which the comparison is done in regards to.
[in]obj2The ProSHADE_data object for which the comparison is done in regards from - the E matrices will be saved into this object.
[in]bandThe bandwidth of the SH value for which this should be done.
[in]orderThe order of the SH value for which this should be done.
[in]obj1ValsAlready allocated array of proshade_double to which integrated values will be saved. It must have size equal to minimum of spheres in the two compared objects.
[in]obj2ValsAlready allocated array of proshade_double to which integrated values will be saved. It must have size equal to minimum of spheres in the two compared objects.
[in]integOrderThe Gauss-Legendre integration order to be used.
[in]abscissasThe pre-computed abscissas for the Gauss-Legendre integration.
[in]weightsThe pre-computed weights for the Gauss-Legendre integration.
[in]integRangeThe range in angstroms between the smalleds and largest shell which are integrated over (might not be 0 to max for progressive shell sampling).
[in]sphereDistThe distance between any two spheres.
[out]sphereRangeThe distance between the smallest and largest usable sphere (usable as in having the required band).

Definition at line 462 of file ProSHADE_distances.cpp.

463 {
464  //================================================ Initialise local values
465  proshade_unsign obj1ValsIter = 0;
466  proshade_unsign obj2ValsIter = 0;
467  proshade_unsign integOrderU = static_cast< proshade_unsign > ( integOrder );
468 
469  //================================================ Set sphere counters
470  proshade_unsign minSphere = std::min( obj1->getMaxSpheres(), obj2->getMaxSpheres() );
471  proshade_unsign maxSphere = 0;
472 
473  //================================================ For each radius, deal with weights
474  for ( proshade_unsign radiusIter = 0; radiusIter < std::min( obj1->getMaxSpheres(), obj2->getMaxSpheres() ); radiusIter++ )
475  {
476  //============================================ Get only values where the shell has the band
477  if ( std::min ( obj1->getShellBandwidth ( radiusIter ), obj2->getShellBandwidth ( radiusIter ) ) <= static_cast< proshade_unsign > ( bandIter ) ) { continue; }
478  minSphere = std::min ( radiusIter, minSphere );
479  maxSphere = std::max ( radiusIter, maxSphere );
480 
481  //============================================ Get the magnitudes for weighting
482  computeSphericalHarmonicsMagnitude ( obj1, bandIter, orderIter, radiusIter, &(obj1Vals[obj1ValsIter]) );
483  computeSphericalHarmonicsMagnitude ( obj2, bandIter, orderIter, radiusIter, &(obj2Vals[obj2ValsIter]) );
484  obj1ValsIter += 1;
485  obj2ValsIter += 1;
486  }
487 
488  //================================================ Integrate weights
489  proshade_single minSphereRad = obj1->getSpherePosValue ( minSphere ) - ( sphereDist * 0.5f );
490  proshade_single maxSphereRad = obj1->getSpherePosValue ( maxSphere ) + ( sphereDist * 0.5f );
491 
492  obj1->setIntegrationWeightCumul ( ProSHADE_internal_maths::gaussLegendreIntegrationReal ( obj1Vals, obj1ValsIter, integOrderU, abscissas, weights, static_cast< proshade_double > ( maxSphereRad - minSphereRad ), static_cast< proshade_double > ( sphereDist ) ) );
493  obj2->setIntegrationWeightCumul ( ProSHADE_internal_maths::gaussLegendreIntegrationReal ( obj2Vals, obj2ValsIter, integOrderU, abscissas, weights, static_cast< proshade_double > ( maxSphereRad - minSphereRad ), static_cast< proshade_double > ( sphereDist ) ) );
494 
495  //================================================ Done
496  return ( static_cast< proshade_double > ( maxSphereRad - minSphereRad ) );
497 
498 }

◆ generateSO3CoeffsFromEMatrices()

void ProSHADE_internal_distances::generateSO3CoeffsFromEMatrices ( ProSHADE_internal_data::ProSHADE_data obj2,
ProSHADE_settings settings 
)

This function converts the E matrices to SO(3) coefficients.

This function starts by allocating the memory for the SO(3) coefficients and their inverse. It then proceeds to convert the E matrix values into the SO(3) transform coefficients by applying the Wigner normalisation factor and changing the sign as required by SOFT library. Upon termination, the coeffs will be saved in the obj2 class.

Parameters
[in]obj2The second ProSHADE_data object which is compared to the first.
[in]settingsA pointer to settings class containing all the information required for the task.

Definition at line 732 of file ProSHADE_distances.cpp.

733 {
734  //================================================ Report progress
735  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, "Converting E matrices to SO(3) coefficients.", settings->messageShift );
736 
737  //================================================ Allocate memory for the coefficients
738  obj2->allocateSO3CoeffsSpace ( obj2->getEMatDim ( ) );
739 
740  //================================================ Initialise local variables
741  proshade_double wigNorm, hlpValReal, hlpValImag;
742  proshade_double signValue = 1.0;
743  proshade_unsign indexO;
744  proshade_complex hlpVal;
745 
746  //================================================ For each band (l)
747  for ( proshade_signed bandIter = 0; bandIter < static_cast<proshade_signed> ( obj2->getEMatDim ( ) ); bandIter++ )
748  {
749  //============================================ Get wigner normalisation factor
750  wigNorm = 2.0 * M_PI * sqrt ( 2.0 / (2.0 * static_cast< proshade_double > ( bandIter ) + 1.0 ) );
751 
752  //============================================ For each order (m)
753  for ( proshade_signed orderIter = 0; orderIter < ( ( bandIter * 2 ) + 1 ); orderIter++ )
754  {
755  //======================================== Set the sign
756  if ( ( orderIter - bandIter + bandIter ) % 2 ) { signValue = -1.0 ; }
757  else { signValue = 1.0 ; }
758 
759  //======================================== For each order2 (m')
760  for ( proshade_signed order2Iter = 0; order2Iter < ( ( bandIter * 2 ) + 1 ); order2Iter++ )
761  {
762  //==================================== Find output index
763  indexO = static_cast< proshade_unsign > ( so3CoefLoc ( static_cast< int > ( orderIter - bandIter ), static_cast< int > ( order2Iter - bandIter ), static_cast< int > ( bandIter ), static_cast< int > ( obj2->getEMatDim ( ) ) ) );
764 
765  //==================================== Compute and save the SO(3) coefficients
766  obj2->getEMatrixValue ( static_cast< proshade_unsign > ( bandIter ), static_cast< proshade_unsign > ( orderIter ), static_cast< proshade_unsign > ( order2Iter ), &hlpValReal, &hlpValImag );
767  hlpVal[0] = hlpValReal * wigNorm * signValue;
768  hlpVal[1] = hlpValImag * wigNorm * signValue;
769  obj2->setSO3CoeffValue ( indexO, hlpVal );
770 
771  //==================================== Switch the sign value
772  signValue *= -1.0;
773  }
774  }
775  }
776 
777  //================================================ Report progress
778  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 3, "SO(3) coefficients obtained.", settings->messageShift );
779 
780  //================================================ Done
781  return ;
782 
783 }

◆ isBandWithinShell()

bool ProSHADE_internal_distances::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.

This function simply checks if a particular sphere bandwidth limit is higher than a requested band, returning true if it is and false if not.

Parameters
[in]bandInQuestionThe value of the band existence of which is to be checked.
[in]shellInQuestionThe index of the shell for which the band existence is checked.
[in]spheresThe ProSHADE structure holding all the shells and their information.

Definition at line 139 of file ProSHADE_distances.cpp.

140 {
141  if ( bandInQuestion < spheres[shellInQuestion]->getLocalBandwidth() )
142  {
143  return ( true );
144  }
145  else
146  {
147  return ( false );
148  }
149 }

◆ normaliseEMatrices()

void ProSHADE_internal_distances::normaliseEMatrices ( ProSHADE_internal_data::ProSHADE_data obj1,
ProSHADE_internal_data::ProSHADE_data obj2,
ProSHADE_settings settings 
)

This function normalises the E matrices.

This function assumes that the E matrices and integration magnitude weights were already computed. It now proceeds to compute the weighting factor (sqrt of the product of the magnitudes of the two objects) and apply it to the E matrices. This normalisation is similar in formula and meaning to the Pearson's correlation coefficient normalisation.

Parameters
[in]obj1The first ProSHADE_data object for which the computation is done.
[in]obj2The second ProSHADE_data object for which the computation is done.
[in]settingsA pointer to settings class containing all the information required for the task.

Definition at line 604 of file ProSHADE_distances.cpp.

605 {
606  //================================================ Report progress
607  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 3, "Starting E matrices normalisation.", settings->messageShift );
608 
609  //================================================ Normalise by the Pearson's c.c. like formula
610  proshade_double eMatNormFactor = std::sqrt ( obj1->getIntegrationWeight() * obj2->getIntegrationWeight() );
611 
612  //================================================ If this is self-correlation (i.e. obj1 == obj2), then divide normalisation factor by 2 as the weight was applied cumulatively!
613  if ( obj1->inputOrder == obj2->inputOrder ) { eMatNormFactor /= 2.0; }
614 
615  for ( proshade_unsign bandIter = 0; bandIter < std::min ( obj1->getMaxBand(), obj2->getMaxBand() ); bandIter++ )
616  {
617  //============================================ For each combination of m and m' for E matrices
618  for ( proshade_unsign orderIter = 0; orderIter < ( ( bandIter * 2 ) + 1 ); orderIter++ )
619  {
620  for ( proshade_unsign order2Iter = 0; order2Iter < ( ( bandIter * 2 ) + 1 ); order2Iter++ )
621  {
622  obj2->normaliseEMatrixValue ( bandIter, orderIter, order2Iter, eMatNormFactor );
623  }
624  }
625  }
626 
627  //================================================ Report progress
628  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 4, "E matrices normalised.", settings->messageShift );
629 
630  //================================================ Done
631  return ;
632 
633 }

◆ prepareInvSOFTPlan()

void ProSHADE_internal_distances::prepareInvSOFTPlan ( fftw_plan *  inverseSO3,
int  band,
fftw_complex *  work1,
proshade_complex *  invCoeffs 
)

This function prepares the FFTW plan for the inverse SO(3) transform.

Parameters
[in]inverseSO3The FFTW_PLAN pointer where the result will be saved.
[in]bandThe bandwidth of the computations.
[in]work1The workspace to be used for the computation.
[in]invCoeffsThe pointer to where the inverse SOFT transform results will be saved.

Definition at line 816 of file ProSHADE_distances.cpp.

817 {
818  //================================================ Prepare the plan describing variables
819  int howmany = 4 * band * band;
820  int idist = 2 * band;
821  int odist = 2 * band;
822  int rank = 2;
823 
824  int inembed[2], onembed[2];
825  inembed[0] = 2 * band;
826  inembed[1] = 4 * band * band;
827  onembed[0] = 2 * band;
828  onembed[1] = 4 * band * band;
829 
830  int istride = 1;
831  int ostride = 1;
832 
833  int na[2];
834  na[0] = 1;
835  na[1] = 2 * band;
836 
837  //================================================ Create the plan
838  *inverseSO3 = fftw_plan_many_dft ( rank,
839  na,
840  howmany,
841  work1,
842  inembed,
843  istride,
844  idist,
845  invCoeffs,
846  onembed,
847  ostride,
848  odist,
849  FFTW_FORWARD,
850  FFTW_ESTIMATE );
851 
852  //================================================ Done
853  return ;
854 
855 }

◆ releaseInvSOFTMemory()

void ProSHADE_internal_distances::releaseInvSOFTMemory ( proshade_complex *&  work1,
proshade_complex *&  work2,
proshade_double *&  work3 
)

This function releases the memory used for computation of the inverse SOFT transform.

Parameters
[in]work1The first workspace pointer to be released.
[in]work2The second workspace pointer to be released.
[in]work3The third workspace pointer to be released.

Definition at line 863 of file ProSHADE_distances.cpp.

864 {
865  //================================================ Release memory
866  delete[] work1;
867  delete[] work2;
868  delete[] work3;
869 
870  //================================================ Done
871  return ;
872 }

◆ releaseTrSigmaWorkspace()

void ProSHADE_internal_distances::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 descriptor computation.

Parameters
[in]obj1ValsArray to hold the shell values for the first object integgration.
[in]obj2ValsArray to hold the shell values for the second object integgration.
[in]GLabscissasAn array to hold the pre-computed anscissas for the Gauss-Legendre integration.
[in]glWeightsAn array to hold the pre-computed weights for the Gauss-Legendre integration.
[in]radiiValsA complex array to hold the results of combining spherical harmonics coefficients of the two objects for each shell.

Definition at line 508 of file ProSHADE_distances.cpp.

509 {
510  //================================================ Release memory
511  delete[] obj1Vals;
512  delete[] obj2Vals;
513  delete[] radiiVals;
514  delete[] GLabscissas;
515  delete[] GLweights;
516 
517  //================================================ Set to NULL
518  obj1Vals = nullptr;
519  obj2Vals = nullptr;
520  radiiVals = nullptr;
521  GLabscissas = nullptr;
522  GLweights = nullptr;
523 
524  //================================================ Done
525  return ;
526 
527 }
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:604
ProSHADE_settings::integOrder
proshade_unsign integOrder
The order required for full Gauss-Legendre integration between the spheres.
Definition: ProSHADE_settings.hpp:70
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:3962
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:714
ProSHADE_settings::computeTraceSigmaDesc
bool computeTraceSigmaDesc
If true, the trace sigma descriptor will be computed, otherwise all its computations will be omitted.
Definition: ProSHADE_settings.hpp:119
ProSHADE_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:624
ProSHADE_settings::computeRotationFuncDesc
bool computeRotationFuncDesc
If true, the rotation function descriptor will be computed, otherwise all its computations will be om...
Definition: ProSHADE_settings.hpp:120
ProSHADE_internal_spheres::ProSHADE_sphere::getLocalBandwidth
proshade_unsign getLocalBandwidth(void)
This function returns the local bandwidth.
Definition: ProSHADE_spheres.cpp:391
ProSHADE_internal_data::ProSHADE_data::sphericalHarmonics
proshade_complex ** sphericalHarmonics
A set of spherical harmonics values arrays for each sphere.
Definition: ProSHADE_data.hpp:121
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:329
ProSHADE_internal_distances::computeEMatricesForLM
void computeEMatricesForLM(ProSHADE_internal_data::ProSHADE_data *obj1, ProSHADE_internal_data::ProSHADE_data *obj2, int bandIter, int orderIter, proshade_complex *radiiVals, int 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:391
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:4305
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:4354
ProSHADE_exception
This class is the representation of ProSHADE exception.
Definition: ProSHADE_exceptions.hpp:37
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:3977
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:816
ProSHADE_internal_data::ProSHADE_data::inputOrder
proshade_unsign inputOrder
This value is the input order - it is useful to know for writing out files, so that they would not ov...
Definition: ProSHADE_data.hpp:140
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:863
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:3727
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:4370
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:3915
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:111
ProSHADE_internal_data::ProSHADE_data::getSO3Coeffs
proshade_complex * getSO3Coeffs(void)
This function allows access to the SO(3) coefficients array.
Definition: ProSHADE_data.cpp:4003
ProSHADE_settings::maxSphereDists
proshade_single maxSphereDists
The distance between spheres in spherical mapping for the largest sphere.
Definition: ProSHADE_settings.hpp:67
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:3926
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:4319
ProSHADE_settings::verbose
proshade_signed verbose
Should the software report on the progress, or just be quiet? Value between -1 (nothing) and 4 (loud)
Definition: ProSHADE_settings.hpp:152
ProSHADE_internal_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:807
ProSHADE_settings::messageShift
proshade_signed messageShift
This value allows shifting the messages to create more readable log for sub-processes.
Definition: ProSHADE_settings.hpp:153
ProSHADE_internal_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:4029
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:3938
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_distances::generateSO3CoeffsFromEMatrices
void generateSO3CoeffsFromEMatrices(ProSHADE_internal_data::ProSHADE_data *obj2, ProSHADE_settings *settings)
This function converts the E matrices to SO(3) coefficients.
Definition: ProSHADE_distances.cpp:732
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_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_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:3774
ProSHADE_internal_data::ProSHADE_data::getMaxBand
proshade_unsign getMaxBand(void)
This function returns the maximum band value for the object.
Definition: ProSHADE_data.cpp:3748
ProSHADE_settings::resolutionOversampling
proshade_single resolutionOversampling
How much (%) should the requested resolution be over-sampled by map re-sampling?
Definition: ProSHADE_settings.hpp:53
ProSHADE_settings::integApproxSteps
proshade_unsign integApproxSteps
The number of steps taken in the approximation of Legendre polynomial decomposition to terms.
Definition: ProSHADE_settings.hpp:71
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:540
ProSHADE_internal_data::ProSHADE_data::getSpherePosValue
proshade_single getSpherePosValue(proshade_unsign shell)
This function allows access to sphere positions.
Definition: ProSHADE_data.cpp:3950
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:3758
ProSHADE_internal_distances::computeInverseSOFTTransform
void computeInverseSOFTTransform(ProSHADE_internal_data::ProSHADE_data *obj2, ProSHADE_settings *settings)
This function computes the inverse SO(3) transform.
Definition: ProSHADE_distances.cpp:885
ProSHADE_internal_distances::computeWeightsForEMatricesForLM
proshade_double computeWeightsForEMatricesForLM(ProSHADE_internal_data::ProSHADE_data *obj1, ProSHADE_internal_data::ProSHADE_data *obj2, int bandIter, int orderIter, proshade_double *obj1Vals, proshade_double *obj2Vals, int 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:462
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:73
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:508
ProSHADE_internal_distances::computeSphericalHarmonicsMagnitude
void computeSphericalHarmonicsMagnitude(ProSHADE_internal_data::ProSHADE_data *obj, int band, int 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:358
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:792
ProSHADE_internal_data::ProSHADE_data::setEMatrixValue
void setEMatrixValue(int band, int order1, int order2, proshade_complex val)
This function allows setting the E matrix value.
Definition: ProSHADE_data.cpp:4336
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:3992
ProSHADE_internal_data::ProSHADE_data::getEMatDim
proshade_unsign getEMatDim(void)
This function allows access to the maximum band for the E matrix.
Definition: ProSHADE_data.cpp:4014
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_data::ProSHADE_data::allocateEMatrices
void allocateEMatrices(proshade_unsign band, proshade_single oversamplingRatio)
This function allocates the required memory for the E matrices.
Definition: ProSHADE_distances.cpp:276
ProSHADE_internal_maths::getLegendreAbscAndWeights
void getLegendreAbscAndWeights(proshade_unsign order, proshade_double *abscissas, proshade_double *weights, proshade_unsign noSteps)
Function to prepare abscissas and weights for Gauss-Legendre integration using the Glaser-Liu-Rokhlin...
Definition: ProSHADE_maths.cpp:289
ProSHADE_settings::noIntegrationSpeedup
bool noIntegrationSpeedup
This option turns off the integration speedup (only using abscissas and weights up to appropriate l f...
Definition: ProSHADE_settings.hpp:72
ProSHADE_internal_data::ProSHADE_data::spheres
ProSHADE_internal_spheres::ProSHADE_sphere ** spheres
The set of concentric spheres to which the intermal density map has been projected.
Definition: ProSHADE_data.hpp:120
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:707
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_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:118
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:117