![]() |
ProSHADE
0.7.6.2 (DEC 2021)
Protein Shape Detection
|
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, 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, weighting it by the radius^2 (for integration). More... | |
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 into the obj2 parameter. More... | |
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 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 *obj1, 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 *obj1, 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... | |
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.
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.
[in] | work1 | The first workspace pointer to be allocated. |
[in] | work2 | The second workspace pointer to be allocated. |
[in] | work3 | The third workspace pointer to be allocated. |
[in] | band | The bandwidth of the computations (this determines how much workspace will be required). |
Definition at line 761 of file ProSHADE_distances.cpp.
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.
[in] | minSpheres | The minima of the number of spheres available in the compared objects. |
[in] | intOrder | The integration order for the computation. |
[in] | obj1Vals | Array to hold the shell values for the first object integgration. |
[in] | obj2Vals | Array to hold the shell values for the second object integgration. |
[in] | GLabscissas | An array to hold the pre-computed anscissas for the Gauss-Legendre integration. |
[in] | glWeights | An array to hold the pre-computed weights for the Gauss-Legendre integration. |
[in] | radiiVals | A complex array to hold the results of combining spherical harmonics coefficients of the two objects for each shell. |
Definition at line 312 of file ProSHADE_distances.cpp.
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!
[in] | obj1 | The first ProSHADE_data object for which the computation is done. |
[in] | obj2 | The second ProSHADE_data object for which the computation is done. |
[in] | settings | A pointer to settings class containing all the information required for the task. |
Definition at line 510 of file ProSHADE_distances.cpp.
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 | ||
) |
This function computes the E matrix un-weighted values for a given band and order and saves these into the obj2 parameter.
[in] | obj1 | The ProSHADE_data object for which the comparison is done in regards to. |
[in] | obj2 | The ProSHADE_data object for which the comparison is done in regards from - the E matrices will be saved into this object. |
[in] | band | The bandwidth of the SH value for which this should be done. |
[in] | order | The order of the SH value for which this should be done. |
[in] | radiiVals | Already 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] | integOrder | The Gauss-Legendre integration order to be used. |
[in] | abscissas | The pre-computed abscissas for the Gauss-Legendre integration. |
[in] | weights | The pre-computed weights for the Gauss-Legendre integration. |
[in] | integRange | The range in angstroms between the smalleds and largest shell which are integrated over (might not be 0 to max for progressive shell sampling). |
[in] | sphereDist | The distance between any two spheres. |
Definition at line 370 of file ProSHADE_distances.cpp.
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.
[in] | obj1 | The first ProSHADE_data object against which comparison is done. |
[in] | obj2 | The second ProSHADE_data object which is compared to the first. |
[in] | settings | A pointer to settings class containing all the information required for the task. |
Definition at line 161 of file ProSHADE_distances.cpp.
void ProSHADE_internal_distances::computeInverseSOFTTransform | ( | ProSHADE_internal_data::ProSHADE_data * | obj1, |
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.
[in] | obj1 | The first ProSHADE_data object against which comparison is done. |
[in] | obj2 | The second ProSHADE_data object which is compared to the first. |
[in] | settings | A pointer to settings class containing all the information required for the task. |
Definition at line 854 of file ProSHADE_distances.cpp.
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.
[in] | obj1 | The first ProSHADE_data object against which comparison is done. |
[in] | obj2 | The second ProSHADE_data object which is compared to the first. |
[in] | settings | A pointer to settings class containing all the information required for the task. |
[out] | ret | The final normalised value of the rotation function descriptor for the two objects. |
Definition at line 908 of file ProSHADE_distances.cpp.
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.
[in] | obj1 | The first ProSHADE_data object against which comparison is done. |
[in] | obj2 | The second ProSHADE_data object which is compared to the first. |
[in] | settings | A pointer to settings class containing all the information required for the task. |
[in] | minCommonBands | The number of common bands between the two objects. |
[in] | minCommonShells | The index of highest common shell in both objects. |
[in] | bandDists | Empty 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.
void ProSHADE_internal_distances::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, weighting it by the radius^2 (for integration).
[in] | obj | The ProSHADE_data object for which the computation is to be done. |
[in] | band | The bandwidth of the SH value for which this should be done. |
[in] | order | The order of the SH value for which this should be done. |
[in] | radius | The shell of the SH value for which this should be done. |
[in] | result | The location where the result is to be saved. |
Definition at line 341 of file ProSHADE_distances.cpp.
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.
[in] | obj1 | The first ProSHADE_data object against which comparison is done. |
[in] | obj2 | The second ProSHADE_data object which is compared to the first. |
[in] | settings | A pointer to settings class containing all the information required for the task. |
[out] | ret | The final normalised value of the trace sigma descriptor for the two objects. |
Definition at line 616 of file ProSHADE_distances.cpp.
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 | ||
) |
This function computes the E matrix weight values for a given band and order and saves these into the appropriate objects.
[in] | obj1 | The ProSHADE_data object for which the comparison is done in regards to. |
[in] | obj2 | The ProSHADE_data object for which the comparison is done in regards from - the E matrices will be saved into this object. |
[in] | band | The bandwidth of the SH value for which this should be done. |
[in] | order | The order of the SH value for which this should be done. |
[in] | obj1Vals | Already 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] | obj2Vals | Already 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] | integOrder | The Gauss-Legendre integration order to be used. |
[in] | abscissas | The pre-computed abscissas for the Gauss-Legendre integration. |
[in] | weights | The pre-computed weights for the Gauss-Legendre integration. |
[in] | integRange | The range in angstroms between the smalleds and largest shell which are integrated over (might not be 0 to max for progressive shell sampling). |
[in] | sphereDist | The distance between any two spheres. |
[out] | sphereRange | The distance between the smallest and largest usable sphere (usable as in having the required band). |
Definition at line 433 of file ProSHADE_distances.cpp.
void ProSHADE_internal_distances::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.
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.
[in] | obj1 | The first ProSHADE_data object against which comparison is done. |
[in] | obj2 | The second ProSHADE_data object which is compared to the first. |
[in] | settings | A pointer to settings class containing all the information required for the task. |
Definition at line 701 of file ProSHADE_distances.cpp.
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.
[in] | bandInQuestion | The value of the band existence of which is to be checked. |
[in] | shellInQuestion | The index of the shell for which the band existence is checked. |
[in] | spheres | The ProSHADE structure holding all the shells and their information. |
Definition at line 139 of file ProSHADE_distances.cpp.
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.
[in] | obj1 | The first ProSHADE_data object for which the computation is done. |
[in] | obj2 | The second ProSHADE_data object for which the computation is done. |
[in] | settings | A pointer to settings class containing all the information required for the task. |
Definition at line 572 of file ProSHADE_distances.cpp.
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.
[in] | inverseSO3 | The FFTW_PLAN pointer where the result will be saved. |
[in] | band | The bandwidth of the computations. |
[in] | work1 | The workspace to be used for the computation. |
[in] | invCoeffs | The pointer to where the inverse SOFT transform results will be saved. |
Definition at line 785 of file ProSHADE_distances.cpp.
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.
[in] | work1 | The first workspace pointer to be released. |
[in] | work2 | The second workspace pointer to be released. |
[in] | work3 | The third workspace pointer to be released. |
Definition at line 832 of file ProSHADE_distances.cpp.
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.
[in] | obj1Vals | Array to hold the shell values for the first object integgration. |
[in] | obj2Vals | Array to hold the shell values for the second object integgration. |
[in] | GLabscissas | An array to hold the pre-computed anscissas for the Gauss-Legendre integration. |
[in] | glWeights | An array to hold the pre-computed weights for the Gauss-Legendre integration. |
[in] | radiiVals | A complex array to hold the results of combining spherical harmonics coefficients of the two objects for each shell. |
Definition at line 478 of file ProSHADE_distances.cpp.