 |
ProSHADE
0.7.6.6 (JUL 2022)
Protein Shape Detection
|
Go to the documentation of this file.
26 #ifndef PROSHADE_SYMMETRY
27 #define PROSHADE_SYMMETRY
37 bool isSymmetrySame ( std::vector< proshade_double* >* ret, proshade_double* sym, proshade_double simThres, proshade_signed* matchedPos );
38 bool isSymmetrySame ( std::vector< proshade_double* >* ret, proshade_double* sym, proshade_double simThres, proshade_signed* matchedPos, proshade_double fscVal );
39 bool detectTetrahedralSymmetry ( std::vector< proshade_double* >* CSymList, proshade_double axErr, proshade_double minPeakHeight );
40 void findTetra4C3s ( std::vector< proshade_double* >* CSymList, std::vector< proshade_double* >* ret,
42 proshade_signed verbose, proshade_signed messageShift, proshade_double minPeakHeight );
44 proshade_double* sym, proshade_double axErr, proshade_double angle,
bool improve,
45 proshade_unsign pos = 0 );
46 bool findMissingAxes ( std::vector< std::vector< proshade_unsign > >* possibilities,
47 std::vector< proshade_double* >* CSymList, proshade_unsign requiredNoAxes,
48 proshade_double axErr, proshade_double angle, proshade_unsign fold,
50 proshade_double
missingAxisHeight ( proshade_double xVal, proshade_double yVal, proshade_double zVal,
52 std::vector < proshade_double* >
findMissingAxisPoints ( proshade_double xVal, proshade_double yVal, proshade_double zVal,
54 bool sortArrVecHlp (
const proshade_double* a,
const proshade_double* b );
55 void saveMissingAxisNewOnly ( std::vector< proshade_double* >* axVec, proshade_double axX, proshade_double axY,
56 proshade_double axZ, proshade_double height, proshade_unsign fold, proshade_double axErr );
58 std::vector< proshade_unsign >* grp, std::vector< proshade_double* >* hlpVec,
59 proshade_double axErr, proshade_double angle, proshade_unsign fold,
60 proshade_double minPeakHeight );
61 void findTetra3C2s ( std::vector< proshade_double* >* CSymList, std::vector< proshade_double* >* ret,
63 proshade_signed verbose, proshade_signed messageShift, proshade_double minPeakHeight );
64 bool testGroupAgainstGroup ( std::vector< proshade_double* >* CSymList, std::vector< proshade_unsign >* grp1,
65 std::vector< proshade_double* >* RetList, std::vector< proshade_unsign >* grp2,
66 proshade_double angle, proshade_double axErr );
67 bool detectOctahedralSymmetry ( std::vector< proshade_double* >* CSymList, proshade_double axErr, proshade_double minPeakHeight );
68 void findOcta3C4s ( std::vector< proshade_double* >* CSymList, std::vector< proshade_double* >* ret,
70 proshade_signed verbose, proshade_signed messageShift, proshade_double minPeakHeight );
71 void findOcta4C3s ( std::vector< proshade_double* >* CSymList, std::vector< proshade_double* >* ret,
73 proshade_signed verbose, proshade_signed messageShift, proshade_double minPeakHeight );
74 void findOcta6C2s ( std::vector< proshade_double* >* CSymList, std::vector< proshade_double* >* ret,
76 proshade_signed verbose, proshade_signed messageShift, proshade_double minPeakHeight );
78 std::vector< proshade_double* >* CSymList, std::vector< proshade_double* >* ret, std::vector<
79 proshade_unsign >* retGroup, proshade_unsign requiredNoAxes, proshade_double axErr,
80 proshade_unsign noMatchesG1, proshade_double angle1, proshade_unsign noMatchesG2,
82 proshade_signed
addAxisUnlessSame ( proshade_unsign fold, proshade_double axX, proshade_double axY, proshade_double axZ,
83 proshade_double axHeight, std::vector< proshade_double* >* prosp, proshade_double axErr );
84 proshade_signed
addAxisUnlessSame ( proshade_unsign fold, proshade_double axX, proshade_double axY, proshade_double axZ,
85 proshade_double axHeight, proshade_double averageFSC, std::vector< proshade_double* >* prosp, proshade_double axErr );
87 proshade_unsign fold, proshade_double axX, proshade_double axY, proshade_double axZ,
88 std::vector< proshade_double* >* prosp, proshade_double axErr, proshade_unsign noMatchesG1,
89 proshade_double angle1, proshade_unsign noMatchesG2, proshade_double angle2,
91 bool detectIcosahedralSymmetry ( std::vector< proshade_double* >* CSymList, proshade_double axErr, proshade_double minPeakHeight );
92 void findIcos6C5s ( std::vector< proshade_double* >* CSymList, std::vector< proshade_double* >* ret,
94 proshade_signed verbose, proshade_signed messageShift, proshade_double minPeakHeight );
95 void findIcos10C3s ( std::vector< proshade_double* >* CSymList, std::vector< proshade_double* >* ret,
97 proshade_signed verbose, proshade_signed messageShift, proshade_double minPeakHeight );
98 void findIcos15C2s ( std::vector< proshade_double* >* CSymList, std::vector< proshade_double* >* ret,
100 proshade_signed verbose, proshade_signed messageShift, proshade_double minPeakHeight );
102 std::vector< proshade_double* >* CSymList, std::vector< proshade_double* >* ret, std::vector<
103 proshade_unsign >* retGroup, proshade_unsign requiredNoAxes, proshade_double axErr,
104 proshade_unsign noMatchesG1, proshade_double angle1, proshade_unsign noMatchesG2,
105 proshade_double angle2, proshade_unsign noMatchesG3, proshade_double angle3,
108 proshade_unsign fold, proshade_double axX, proshade_double axY, proshade_double axZ,
109 std::vector< proshade_double* >* prosp, proshade_double axErr,
110 proshade_unsign noMatchesG1, proshade_double angle1, proshade_unsign noMatchesG2,
111 proshade_double angle2, proshade_unsign noMatchesG3, proshade_double angle3,
118 void predictIcosAxes ( std::vector< proshade_double* >* CSymList, std::vector< std::vector< proshade_double* > >* ret, proshade_double axErr, proshade_double minPeakHeight );
119 void predictOctaAxes ( std::vector< proshade_double* >* CSymList, std::vector< std::vector< proshade_double* > >* ret, proshade_double axErr, proshade_double minPeakHeight );
120 void predictTetraAxes ( std::vector< proshade_double* >* CSymList, std::vector< std::vector< proshade_double* > >* ret, proshade_double axErr, proshade_double minPeakHeight );
121 std::vector< proshade_unsign >
findReliableUnphasedSymmetries ( std::vector < proshade_double* >* allCs, std::vector < std::vector < proshade_double* > >* allDs, proshade_signed verbose, proshade_signed messageShift, proshade_double tolerance );
122 void allocateCentreOfMapFourierTransforms ( proshade_unsign xDim, proshade_unsign yDim, proshade_unsign zDim, fftw_complex *&origMap, fftw_complex *&origCoeffs, fftw_complex *&rotMapComplex,
123 fftw_complex *&rotCoeffs, fftw_complex *&trFunc, fftw_complex *&trFuncCoeffs, fftw_plan *planForwardFourier, fftw_plan *planForwardFourierRot,
124 fftw_plan *planReverseFourierComb );
125 void releaseCentreOfMapFourierTransforms ( fftw_complex *origMap, fftw_complex *origCoeffs, fftw_complex *rotMapComplex, fftw_complex *rotCoeffs, fftw_complex *trFunc, fftw_complex *trFuncCoeffs,
126 fftw_plan planForwardFourier, fftw_plan planForwardFourierRot, fftw_plan planReverseFourierComb );
128 fftw_complex* rotMapComplex, fftw_complex* rotCoeffs, fftw_plan planForwardFourierRot, fftw_complex* trFuncCoeffs,
129 fftw_complex* trFunc, fftw_plan planReverseFourierComb );
131 fftw_complex* rotMapComplex, fftw_complex* rotCoeffs, fftw_plan planForwardFourierRot, fftw_complex* trFuncCoeffs,
132 fftw_complex* trFunc, fftw_plan planReverseFourierComb );
void findTetra3C2s(std::vector< proshade_double * > *CSymList, std::vector< proshade_double * > *ret, proshade_double axErr, ProSHADE_internal_data::ProSHADE_data *dataObj, proshade_signed verbose, proshade_signed messageShift, proshade_double minPeakHeight)
This function takes the list of C symmetries and finds the 3 C2 symmetries with correct angles requir...
bool detectOctahedralSymmetry(std::vector< proshade_double * > *CSymList, proshade_double axErr, proshade_double minPeakHeight)
This function takes the list of C symmetries and decides whether basic requirements for octahhedral s...
This is the header file containing declarations of functions required for computation of shape distan...
bool findMissingAxes(std::vector< std::vector< proshade_unsign > > *possibilities, std::vector< proshade_double * > *CSymList, proshade_unsign requiredNoAxes, proshade_double axErr, proshade_double angle, proshade_unsign fold, ProSHADE_internal_data::ProSHADE_data *dataObj, proshade_double minPeakHeight)
This function tries to find an axis which would complete a particular group of axes for polyhedral sy...
bool sortArrVecHlp(const proshade_double *a, const proshade_double *b)
This function compares two arrays of two based on the first number.
void findIcos6C5s(std::vector< proshade_double * > *CSymList, std::vector< proshade_double * > *ret, proshade_double axErr, ProSHADE_internal_data::ProSHADE_data *dataObj, proshade_signed verbose, proshade_signed messageShift, proshade_double minPeakHeight)
This function takes the list of C symmetries and finds the six C5 symmetries with given angles requir...
bool testGroupAgainstGroup(std::vector< proshade_double * > *CSymList, std::vector< proshade_unsign > *grp1, std::vector< proshade_double * > *RetList, std::vector< proshade_unsign > *grp2, proshade_double angle, proshade_double axErr)
This function compares two groups of axes for a single pair having the required angle.
proshade_signed addAxisUnlessSame(proshade_unsign fold, proshade_double axX, proshade_double axY, proshade_double axZ, proshade_double axHeight, proshade_double averageFSC, std::vector< proshade_double * > *prosp, proshade_double axErr)
This function simply creates a new axis from information in aruments and tests if no such axis alread...
void findPredictedAxesHeights(std::vector< proshade_double * > *ret, ProSHADE_internal_data::ProSHADE_data *dataObj, ProSHADE_settings *settings)
This function finds the rotation function value for all axes supplied in the ret parameter.
std::vector< proshade_double > findPointFromTranslations(ProSHADE_internal_data::ProSHADE_data *symStr, std::vector< std::vector< proshade_double > > symElems, fftw_complex *origCoeffs, fftw_complex *rotMapComplex, fftw_complex *rotCoeffs, fftw_plan planForwardFourierRot, fftw_complex *trFuncCoeffs, fftw_complex *trFunc, fftw_plan planReverseFourierComb)
This function computes the average of optimal translations for a cyclic point group.
void allocateCentreOfMapFourierTransforms(proshade_unsign xDim, proshade_unsign yDim, proshade_unsign zDim, fftw_complex *&origMap, fftw_complex *&origCoeffs, fftw_complex *&rotMapComplex, fftw_complex *&rotCoeffs, fftw_complex *&trFunc, fftw_complex *&trFuncCoeffs, fftw_plan *planForwardFourier, fftw_plan *planForwardFourierRot, fftw_plan *planReverseFourierComb)
This function allocates the required memory for the Fourier transforms required to find the centre of...
This class contains all inputed and derived data for a single structure.
proshade_double missingAxisHeight(proshade_double xVal, proshade_double yVal, proshade_double zVal, ProSHADE_internal_data::ProSHADE_data *dataObj, proshade_unsign fold, proshade_double axErr)
This function searches for the highest peaks average that would produce the required axis and fold.
void searchMissingSymmetrySpace(ProSHADE_internal_data::ProSHADE_data *dataObj, std::vector< proshade_double * > *CSymList, std::vector< proshade_unsign > *grp, std::vector< proshade_double * > *hlpVec, proshade_double axErr, proshade_double angle, proshade_unsign fold, proshade_double minPeakHeight)
This function tests feasible axes against the missing axis criteria, returning a set of matching axes...
bool findMissingAxesDual(std::vector< proshade_unsign > *possibilities, std::vector< proshade_double * > *CSymList, std::vector< proshade_double * > *ret, std::vector< proshade_unsign > *retGroup, proshade_unsign requiredNoAxes, proshade_double axErr, proshade_unsign noMatchesG1, proshade_double angle1, proshade_unsign noMatchesG2, proshade_double angle2, proshade_unsign fold, ProSHADE_internal_data::ProSHADE_data *dataObj)
This function tries to find a particular symmetry axes which would complete a group of symmetries wit...
void optimiseDGroupAngleFromAxesHeights(std::vector< std::vector< proshade_double > > *ret, ProSHADE_internal_data::ProSHADE_data *dataObj, ProSHADE_settings *settings)
This function takes two axes with almost dihedral angle and optimises their relative positions as wel...
bool isSymmetrySame(std::vector< proshade_double * > *ret, proshade_double *sym, proshade_double simThres, proshade_signed *matchedPos)
This function checks if a very similar symmetry is not already saved.
void saveMissingAxisNewOnly(std::vector< proshade_double * > *axVec, proshade_double axX, proshade_double axY, proshade_double axZ, proshade_double height, proshade_unsign fold, proshade_double axErr)
This function saves the recovered information about missing axis into a full symmetry,...
void findIcos10C3s(std::vector< proshade_double * > *CSymList, std::vector< proshade_double * > *ret, proshade_double axErr, ProSHADE_internal_data::ProSHADE_data *dataObj, proshade_signed verbose, proshade_signed messageShift, proshade_double minPeakHeight)
This function takes the list of C symmetries and finds the ten C3 symmetries with correct angles requ...
This namespace contains the symmetry detection related code.
bool detectIcosahedralSymmetry(std::vector< proshade_double * > *CSymList, proshade_double axErr, proshade_double minPeakHeight)
This function takes the list of C symmetries and decides whether basic requirements for isosahedral s...
This class stores all the settings and is passed to the executive classes instead of a multitude of p...
void findOcta4C3s(std::vector< proshade_double * > *CSymList, std::vector< proshade_double * > *ret, proshade_double axErr, ProSHADE_internal_data::ProSHADE_data *dataObj, proshade_signed verbose, proshade_signed messageShift, proshade_double minPeakHeight)
This function takes the list of C symmetries and finds the four C3 symmetries with correct angles req...
std::vector< proshade_unsign > findReliableUnphasedSymmetries(std::vector< proshade_double * > *allCs, std::vector< std::vector< proshade_double * > > *allDs, proshade_signed verbose, proshade_signed messageShift, proshade_double tolerance)
This function checks the list of detected axes (presumably from phaseless symmetry detection) and ret...
std::vector< proshade_double > findTranslationBetweenRotatedAndOriginalMap(ProSHADE_internal_data::ProSHADE_data *symStr, std::vector< proshade_double > symElem, fftw_complex *origCoeffs, fftw_complex *rotMapComplex, fftw_complex *rotCoeffs, fftw_plan planForwardFourierRot, fftw_complex *trFuncCoeffs, fftw_complex *trFunc, fftw_plan planReverseFourierComb)
This function takes a single rotation matrix and procceds to compute the optimal translation between ...
void releaseCentreOfMapFourierTransforms(fftw_complex *origMap, fftw_complex *origCoeffs, fftw_complex *rotMapComplex, fftw_complex *rotCoeffs, fftw_complex *trFunc, fftw_complex *trFuncCoeffs, fftw_plan planForwardFourier, fftw_plan planForwardFourierRot, fftw_plan planReverseFourierComb)
This function releases the allocated memory for the Fourier transforms used to find the centre of the...
bool findMissingAxesTriple(std::vector< proshade_unsign > *possibilities, std::vector< proshade_double * > *CSymList, std::vector< proshade_double * > *ret, std::vector< proshade_unsign > *retGroup, proshade_unsign requiredNoAxes, proshade_double axErr, proshade_unsign noMatchesG1, proshade_double angle1, proshade_unsign noMatchesG2, proshade_double angle2, proshade_unsign noMatchesG3, proshade_double angle3, proshade_unsign fold, ProSHADE_internal_data::ProSHADE_data *dataObj)
This function tries to find a particular symmetry axis which would complete a group of symmetries wit...
void predictOctaAxes(std::vector< proshade_double * > *CSymList, std::vector< std::vector< proshade_double * > > *ret, proshade_double axErr, proshade_double minPeakHeight)
This function predicts all octahedral point group symmetry axes from the cyclic point groups list.
void findOcta6C2s(std::vector< proshade_double * > *CSymList, std::vector< proshade_double * > *ret, proshade_double axErr, ProSHADE_internal_data::ProSHADE_data *dataObj, proshade_signed verbose, proshade_signed messageShift, proshade_double minPeakHeight)
This function takes the list of C symmetries and finds the six C2 symmetries with correct angles requ...
bool testGroupAgainstSymmetry(std::vector< proshade_double * > *CSymList, std::vector< proshade_unsign > *grp, proshade_double *sym, proshade_double axErr, proshade_double angle, bool improve, proshade_unsign pos=0)
This function tests whether a symmetry has particular angle to all members of a group.
void findOcta3C4s(std::vector< proshade_double * > *CSymList, std::vector< proshade_double * > *ret, proshade_double axErr, ProSHADE_internal_data::ProSHADE_data *dataObj, proshade_signed verbose, proshade_signed messageShift, proshade_double minPeakHeight)
This function takes the list of C symmetries and finds the 3 C4 symmetries with perpendicular angles ...
void predictTetraAxes(std::vector< proshade_double * > *CSymList, std::vector< std::vector< proshade_double * > > *ret, proshade_double axErr, proshade_double minPeakHeight)
This function predicts all tetrahedral point group symmetry axes from the cyclic point groups list.
void findTetra4C3s(std::vector< proshade_double * > *CSymList, std::vector< proshade_double * > *ret, proshade_double axErr, ProSHADE_internal_data::ProSHADE_data *dataObj, proshade_signed verbose, proshade_signed messageShift, proshade_double minPeakHeight)
This function takes the list of C symmetries and finds the 4 C3 symmetries with correct angles requir...
proshade_double findPredictedSingleAxisHeight(proshade_double *axis, proshade_double fold, ProSHADE_internal_data::ProSHADE_data *dataObj, ProSHADE_settings *settings)
This function finds the rotation function value for a single axis.
void checkFittingAxisTripleAndSave(std::vector< proshade_unsign > *retGroup, std::vector< proshade_double * > *ret, proshade_unsign fold, proshade_double axX, proshade_double axY, proshade_double axZ, std::vector< proshade_double * > *prosp, proshade_double axErr, proshade_unsign noMatchesG1, proshade_double angle1, proshade_unsign noMatchesG2, proshade_double angle2, proshade_unsign noMatchesG3, proshade_double angle3, ProSHADE_internal_data::ProSHADE_data *dataObj)
This function takes a newly detected "missing" axis and tests it for belonging to the group,...
void findIcos15C2s(std::vector< proshade_double * > *CSymList, std::vector< proshade_double * > *ret, proshade_double axErr, ProSHADE_internal_data::ProSHADE_data *dataObj, proshade_signed verbose, proshade_signed messageShift, proshade_double minPeakHeight)
This function takes the list of C symmetries and finds the fifteen C3 symmetries with correct angles ...
void predictIcosAxes(std::vector< proshade_double * > *CSymList, std::vector< std::vector< proshade_double * > > *ret, proshade_double axErr, proshade_double minPeakHeight)
This function predicts all possible icosahedral point groups symmetry axes from the cyclic point grou...
bool detectTetrahedralSymmetry(std::vector< proshade_double * > *CSymList, proshade_double axErr, proshade_double minPeakHeight)
This function takes the list of C symmetries and decides whether basic requirements for tetrahedral s...
bool checkFittingAxisDualAndSave(std::vector< proshade_unsign > *retGroup, std::vector< proshade_double * > *ret, proshade_unsign fold, proshade_double axX, proshade_double axY, proshade_double axZ, std::vector< proshade_double * > *prosp, proshade_double axErr, proshade_unsign noMatchesG1, proshade_double angle1, proshade_unsign noMatchesG2, proshade_double angle2, ProSHADE_internal_data::ProSHADE_data *dataObj)
This function takes a newly detected "missing" axis and tests it for belonging to the group,...
std::vector< proshade_double * > findMissingAxisPoints(proshade_double xVal, proshade_double yVal, proshade_double zVal, ProSHADE_internal_data::ProSHADE_data *dataObj, proshade_double axErr)
This function searches for all the self-rotation map points conforming to the axis,...