 |
ProSHADE
0.7.6.6 (JUL 2022)
Protein Shape Detection
|
Go to the documentation of this file.
136 std::vector<ProSHADE_internal_spheres::ProSHADE_rotFun_sphere*> sphereMappedRotFun;
156 void readInMAP (
ProSHADE_settings* settings, proshade_double* maskArr =
nullptr, proshade_unsign maskXDim = 0, proshade_unsign maskYDim = 0,
157 proshade_unsign maskZDim = 0, proshade_double* weightsArr =
nullptr, proshade_unsign weigXDim = 0, proshade_unsign weigYDim = 0,
158 proshade_unsign weigZDim = 0 );
166 ProSHADE_data ( std::string strName,
double *mapVals,
int len, proshade_single xDmSz, proshade_single yDmSz,
167 proshade_single zDmSz, proshade_unsign xDmInd, proshade_unsign yDmInd, proshade_unsign zDmInd, proshade_signed xFr,
168 proshade_signed yFr, proshade_signed zFr, proshade_signed xT, proshade_signed yT, proshade_signed zT,
169 proshade_unsign inputO );
173 void readInStructure ( std::string fName, proshade_unsign inputO,
ProSHADE_settings* settings, proshade_double* maskArr =
nullptr, proshade_unsign maskXDim = 0,
174 proshade_unsign maskYDim = 0, proshade_unsign maskZDim = 0, proshade_double* weightsArr =
nullptr, proshade_unsign weigXDim = 0,
175 proshade_unsign weigYDim = 0, proshade_unsign weigZDim = 0 );
177 void writeMap ( std::string fName, std::string title =
"Created by ProSHADE and written by GEMMI",
int mode = 2 );
178 void writePdb ( std::string fName, proshade_double euA = 0.0, proshade_double euB = 0.0, proshade_double euG = 0.0,
179 proshade_double trsX = 0.0, proshade_double trsY = 0.0, proshade_double trsZ = 0.0, proshade_double rotX = 0.0,
180 proshade_double rotY = 0.0, proshade_double rotZ = 0.0,
bool firstModel =
true );
181 void writeGemmi ( std::string fName, gemmi::Structure gemmiStruct, proshade_double euA = 0.0, proshade_double euB = 0.0, proshade_double euG = 0.0,
182 proshade_double trsX = 0.0, proshade_double trsY = 0.0, proshade_double trsZ = 0.0, proshade_double rotX = 0.0,
183 proshade_double rotY = 0.0, proshade_double rotZ = 0.0,
bool firstModel =
true );
184 void writeMask ( std::string fName, proshade_double* mask );
206 bool shellBandExists ( proshade_unsign shell, proshade_unsign bandVal );
208 void allocateEMatrices ( proshade_unsign band, proshade_single oversamplingRatio );
224 int so3CoeffsArrayIndex ( proshade_signed order1, proshade_signed order2, proshade_signed band );
227 std::vector< proshade_double* >* CSyms, proshade_double tolerance, proshade_signed*& cutIndices, fftw_complex*& fCoeffsCut,
228 proshade_signed noBins, proshade_double**& bindata, proshade_signed*& binCounts, proshade_double*& fscByBin,
229 proshade_signed xDim, proshade_signed yDim, proshade_signed zDim );
231 proshade_signed*& binCounts, proshade_double*& fscByBin, proshade_signed xDim, proshade_signed yDim, proshade_signed zDim );
233 proshade_signed*& binCounts, proshade_double*& fscByBin, proshade_signed xDim, proshade_signed yDim, proshade_signed zDim );
235 proshade_signed*& binCounts, proshade_double*& fscByBin, proshade_signed xDim, proshade_signed yDim, proshade_signed zDim );
245 void prepareFSCFourierMemory ( proshade_signed*& cutIndices, fftw_complex*& fCoeffsCut, proshade_signed* noBins, proshade_double**& bindata,
246 proshade_signed*& binCounts, proshade_double*& fscByBin, proshade_single resolution,
247 proshade_signed* cutXDim, proshade_signed* cutYDim, proshade_signed* cutZDim );
248 proshade_double
computeFSC (
ProSHADE_settings* settings, std::vector< proshade_double* >* CSym,
size_t symIndex, proshade_signed*& cutIndices,
249 fftw_complex*& fCoeffsCut, proshade_signed noBins, proshade_double**& bindata, proshade_signed*& binCounts,
250 proshade_double*& fscByBin, proshade_signed xDim, proshade_signed yDim, proshade_signed zDim, proshade_unsign rotNumber = 1 );
251 proshade_double
computeFSC (
ProSHADE_settings* settings, proshade_double* sym, proshade_signed*& cutIndices, fftw_complex*& fCoeffsCut, proshade_signed noBins,
252 proshade_double**& bindata, proshade_signed*& binCounts, proshade_double*& fscByBin,
253 proshade_signed xDim, proshade_signed yDim, proshade_signed zDim, proshade_unsign rotNumber = 1 );
255 proshade_signed*& binCounts, proshade_double*& fscByBin, proshade_signed xDim, proshade_signed yDim, proshade_signed zDim );
257 proshade_double**& bindata, proshade_signed*& binCounts, proshade_double*& fscByBin, proshade_signed cutXDim, proshade_signed cutYDim, proshade_signed cutZDim );
258 std::vector<std::vector< proshade_double > >
getAllGroupElements ( std::vector< proshade_unsign > axesList, std::string groupType =
"", proshade_double matrixTolerance = 0.05 );
259 std::vector<std::vector< proshade_double > >
getAllGroupElements ( std::vector < std::vector< proshade_double > >* allCs, std::vector< proshade_unsign > axesList, std::string groupType =
"", proshade_double matrixTolerance = 0.05 );
262 proshade_double**& bindata, proshade_signed*& binCounts, proshade_double*& fscByBin, proshade_signed cutXDim,
263 proshade_signed cutYDim, proshade_signed cutZDim );
270 void zeroPaddToDims ( proshade_unsign xDim, proshade_unsign yDim, proshade_unsign zDim );
272 std::vector< proshade_double >
rotateMapRealSpace ( proshade_double axX, proshade_double axY, proshade_double axZ, proshade_double axAng, proshade_double*& map );
273 std::vector< proshade_double >
rotateMapRealSpaceInPlace ( proshade_double eulA, proshade_double eulB, proshade_double eulG );
274 void rotateFourierCoeffs ( proshade_double axX, proshade_double axY, proshade_double axZ, proshade_double axAng, fftw_complex*& coeffs, fftw_complex*& rotCoeffs,
275 proshade_signed xDim, proshade_signed yDim, proshade_signed zDim );
276 void translateMap ( proshade_double trsX, proshade_double trsY, proshade_double trsZ );
284 std::vector< proshade_double >* ultimateTranslation );
286 std::vector < proshade_double >* finalTranslation );
289 void deepCopyMap ( proshade_double*& saveTo, proshade_signed verbose );
292 proshade_double
getMapValue ( proshade_unsign pos );
295 proshade_double*
getRealSphHarmValue ( proshade_unsign band, proshade_unsign order, proshade_unsign shell );
296 proshade_double*
getImagSphHarmValue ( proshade_unsign band, proshade_unsign order, proshade_unsign shell );
297 proshade_double
getRRPValue ( proshade_unsign band, proshade_unsign sh1, proshade_unsign sh2 );
303 void getEMatrixValue ( proshade_unsign band, proshade_unsign order1, proshade_unsign order2, proshade_double* valueReal, proshade_double* valueImag );
307 void getWignerMatrixValue ( proshade_unsign band, proshade_unsign order1, proshade_unsign order2, proshade_double* valueReal, proshade_double* valueImag );
311 proshade_unsign
getXDim (
void );
312 proshade_unsign
getYDim (
void );
313 proshade_unsign
getZDim (
void );
328 std::vector< std::vector< proshade_double* > >*
getDihedralAxes (
void );
334 void setEMatrixValue (
int band,
int order1,
int order2, proshade_complex val );
335 void normaliseEMatrixValue ( proshade_unsign band, proshade_unsign order1, proshade_unsign order2, proshade_double normF );
337 void setWignerMatrixValue ( proshade_complex val, proshade_unsign band, proshade_unsign order1, proshade_unsign order2 );
341 std::vector<std::vector< proshade_double > >
computeGroupElementsForGroup ( proshade_double xAx, proshade_double yAx, proshade_double zAx, proshade_signed fold );
344 std::vector<std::vector< proshade_double > >* second,
345 proshade_double matrixTolerance,
proshade_double originalMapYCom
The COM of the first map to be loaded/computed without any furhter changes being reflacted along the ...
void computeSphericalHarmonics(ProSHADE_settings *settings)
This function computes the spherical harmonics decomposition for the whole structure.
proshade_double mapMovFromsChangeX
When the map is translated, the xFrom and xTo values are changed. This variable holds how much they h...
proshade_complex ** getEMatrixByBand(proshade_unsign band)
This function allows access to E matrix for a particular band.
void allocateRotatedSHMemory(void)
This function allocates the memory required for storing the rotated Spherical Harmonics coefficients.
std::vector< proshade_double * > getCyclicAxesCopy(void)
This function allows access to the list of detected cyclic axes in non-overwrite fashion.
std::vector< proshade_double * > cyclicSymmetries
This is where the detected cyclic ("C") symmetries will be kept.
void getOverlayRotationFunction(ProSHADE_settings *settings, ProSHADE_internal_data::ProSHADE_data *obj2)
This function computes the overlay rotation function (i.e. the correlation function in SO(3) space).
std::string fileName
This is the original file from which the data were obtained.
proshade_signed zFrom
This is the starting index along the z axis.
proshade_unsign xDimIndices
This is the size of the map cell x dimension in indices.
proshade_single getZDimSize(void)
This function allows access to the map size in angstroms along the Z axis.
std::vector< proshade_double > getBestTranslationMapPeaksAngstrom(ProSHADE_internal_data::ProSHADE_data *staticStructure)
This function gets the optimal translation vector and returns it as a standard library vector....
void readInMAP(ProSHADE_settings *settings, proshade_double *maskArr=nullptr, proshade_unsign maskXDim=0, proshade_unsign maskYDim=0, proshade_unsign maskZDim=0, proshade_double *weightsArr=nullptr, proshade_unsign weigXDim=0, proshade_unsign weigYDim=0, proshade_unsign weigZDim=0)
Function for reading map data using gemmi library.
proshade_single zDimSizeOriginal
This is the size of the map cell z dimension in Angstroms.
void determineRecommendedSymmetry(ProSHADE_settings *settings, proshade_double threshold, proshade_signed *&cutIndices, fftw_complex *&fCoeffsCut, proshade_signed noBins, proshade_double **&bindata, proshade_signed *&binCounts, proshade_double *&fscByBin, proshade_signed xDim, proshade_signed yDim, proshade_signed zDim)
This function takes all the detected symmetry results and decides on which are to be recommended for ...
void createNewMapFromBounds(ProSHADE_settings *settings, ProSHADE_data *&newStr, proshade_signed *newBounds)
This function creates a new structure from the calling structure and new bounds values.
proshade_unsign xGridIndices
As far as I know, this is identical to the xDimIndices.
void computeRotationFunction(ProSHADE_settings *settings)
This function computes the self-rotation function for this structure.
proshade_complex * translationMap
This is where the translation map will be held, if at all used.
std::vector< std::vector< proshade_double > > joinElementsFromDifferentGroups(std::vector< std::vector< proshade_double > > *first, std::vector< std::vector< proshade_double > > *second, proshade_double matrixTolerance, bool combine)
This function joins two group element lists using only unique elements.
proshade_unsign xAxisOrder
This is the order of the x axis.
std::vector< proshade_double * > getRecommendedSymmetryValues(void)
This function simply returns the detected recommended symmetry axes.
void getPredictedIcosahedralSymmetriesList(ProSHADE_settings *settings, std::vector< proshade_double * > *CSymList, proshade_signed *&cutIndices, fftw_complex *&fCoeffsCut, proshade_signed noBins, proshade_double **&bindata, proshade_signed *&binCounts, proshade_double *&fscByBin, proshade_signed xDim, proshade_signed yDim, proshade_signed zDim)
This function predicts a list of all I symmetry axes from the already computed C symmetries list.
proshade_complex ** rotSphericalHarmonics
A set of rotated spherical harmonics values arrays for each sphere, used only if map rotation is requ...
void getReBoxBoundaries(ProSHADE_settings *settings, proshade_signed *&ret)
This function finds the boundaries enclosing positive map values and adds some extra space.
std::vector< proshade_double * > decidePolyFromList(ProSHADE_settings *settings, std::vector< std::vector< proshade_double * > > *polyList, size_t fullGroupSize, std::vector< proshade_double * > *CSyms, proshade_double tolerance, proshade_signed *&cutIndices, fftw_complex *&fCoeffsCut, proshade_signed noBins, proshade_double **&bindata, proshade_signed *&binCounts, proshade_double *&fscByBin, proshade_signed xDim, proshade_signed yDim, proshade_signed zDim)
This function takes a list of predicted polyheral groups and decides which is most likely using the F...
proshade_complex ** sphericalHarmonics
A set of spherical harmonics values arrays for each sphere.
void setIntegrationWeight(proshade_double intW)
This function allows setting the integration weight for the object.
void interpolateMapFromSpheres(proshade_double *&densityMapRotated)
This function interpolates the density map from the sphere mapped data.
proshade_single cAngle
This is the angle c of the map cell in degrees.
std::vector< std::vector< proshade_double * > > * getDihedralAxes(void)
This function allows access to the list of detected dihedral axes.
proshade_complex * so3CoeffsInverse
The inverse coefficients obtained by inverse SO(3) Fourier Transform (SOFT) - i.e....
proshade_signed * getXAxisOrigin(void)
This function allows access to the map X axis origin value.
void normaliseEMatrixValue(proshade_unsign band, proshade_unsign order1, proshade_unsign order2, proshade_double normF)
This function allows normalising the E matrix value.
proshade_double yCom
The COM of the map after processing along the Y-axis.
void getImagTranslationFunction(double *trsFunImag, int len)
This function fills the input array with the imaginary translation function values.
void setWignerMatrixValue(proshade_complex val, proshade_unsign band, proshade_unsign order1, proshade_unsign order2)
This function allows setting the Wigner D matrix value by its band, order1 and order2 co-ordinate.
void allocateRRPMemory()
This function allocates the required memory for the RRP matrices.
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.
proshade_single getXDimSize(void)
This function allows access to the map size in angstroms along the X axis.
proshade_double *& getInternalMap(void)
This function allows access to the first map array value address.
void allocateWignerMatricesSpace()
This function allocates the memory for the Wigner matrices for the calling object.
This class contains all inputed and derived data for a single structure.
proshade_unsign noSpheres
The number of spheres with map projected onto them.
proshade_single getYDimSize(void)
This function allows access to the map size in angstroms along the Y axis.
proshade_unsign inputOrder
This value is the input order - it is useful to know for writing out files, so that they would not ov...
proshade_single zDimSize
This is the size of the map cell z dimension in Angstroms.
proshade_complex * so3Coeffs
The coefficients obtained by SO(3) Fourier Transform (SOFT), in this case derived from the E matrices...
proshade_unsign getMaxSpheres(void)
This function returns the number of spheres which contain the whole object.
void setSO3CoeffValue(proshade_unsign position, proshade_complex val)
This function allows setting the SOFT coefficient values using array position and value.
proshade_double mapCOMProcessChangeX
The change in X axis between the creation of the structure (originalMapXCom) and just before rotation...
proshade_double originalMapXCom
The COM of the first map to be loaded/computed without any furhter changes being reflacted along the ...
void getRotMatrixFromRotFunInds(proshade_signed aI, proshade_signed bI, proshade_signed gI, double *rotMat, int len)
This function takes rotation function indices, converts them to Euler angles and these to rotation ma...
std::vector< proshade_double > getBestRotationMapPeaksEulerAngles(ProSHADE_settings *settings)
This function returns a vector of three floats, the three Euler angles of the best peak in the rotati...
void computeRRPMatrices(ProSHADE_settings *settings)
This function pre-computes the RRP matrices for a data object.
proshade_unsign xDimIndicesOriginal
This is the size of the map cell x dimension in indices.
proshade_double getAnySphereRadius(proshade_unsign shell)
This function allows access to the radius of any particular sphere.
proshade_complex * getSO3Coeffs(void)
This function allows access to the SO(3) coefficients array.
void saveRequestedSymmetryD(ProSHADE_settings *settings, proshade_signed *&cutIndices, fftw_complex *&fCoeffsCut, proshade_signed noBins, proshade_double **&bindata, proshade_signed *&binCounts, proshade_double *&fscByBin, proshade_signed cutXDim, proshade_signed cutYDim, proshade_signed cutZDim)
This function takes the D symmetries and searched for the requested symmetry.
proshade_unsign maxShellBand
The maximum band for any shell of the object.
proshade_double originalPdbRotCenY
The centre of rotation as it relates to the original PDB positions (and not the ProSHADE internal map...
void getSpherePositions(ProSHADE_settings *settings)
This function determines the sphere positions (radii) for sphere mapping.
void invertMirrorMap(ProSHADE_settings *settings)
Function for inverting the map to its mirror image.
std::vector< proshade_double > getMapCOMProcessChange(void)
This function allows access to the translation caused by structure processing.
proshade_double mapCOMProcessChangeZ
The change in Z axis between the creation of the structure (originalMapZCom) and just before rotation...
proshade_single xDimSizeOriginal
This is the size of the map cell x dimension in Angstroms.
proshade_unsign getRecommendedSymmetryFold(void)
This function simply returns the detected recommended symmetry fold.
proshade_double getIntegrationWeight(void)
This function allows access to the integration weight for the object.
void setIntegrationWeightCumul(proshade_double intW)
This function allows setting the cumulative integration weight for the object.
proshade_signed zAxisOriginOriginal
This is the origin position along the z axis.
std::vector< std::vector< proshade_double * > > getDihedralAxesCopy(void)
This function allows access to the list of detected dihedral axes in a non-over-write fashion.
std::vector< proshade_double * > icosahedralSymmetries
This is where the detected icosahedral ("I") symmetries will be kept.
This class contains all inputed and derived data for a single sphere.
std::vector< std::string > getSymmetryAxis(ProSHADE_settings *settings, proshade_unsign axisNo)
This function returns a single symmetry axis as a vector of strings from the recommended symmetry axe...
proshade_complex *** eMatrices
The trace sigma and full rotation function c*conj(c) integral tables.
void centreMapOnCOM(ProSHADE_settings *settings)
This function shits the map so that its COM is in the centre of the map.
proshade_signed * getXFromPtr(void)
This function allows access to the map start along the X axis.
std::vector< std::vector< proshade_double > > getAllGroupElements(std::vector< proshade_unsign > axesList, std::string groupType="", proshade_double matrixTolerance=0.05)
This function returns the group elements as rotation matrices of any defined point group.
proshade_unsign yDimIndices
This is the size of the map cell y dimension in indices.
proshade_signed xAxisOriginOriginal
This is the origin position along the x axis.
void deepCopyMap(proshade_double *&saveTo, proshade_signed verbose)
This function copies the internal map into the supplied pointer, which it also allocates.
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.
void readInGemmi(gemmi::Structure *gemmiStruct, ProSHADE_settings *settings)
Function for reading co-ordinate data from Gemmi object.
proshade_unsign getShellBandwidth(proshade_unsign shell)
This function allows access to the bandwidth of a particular shell.
proshade_unsign yGridIndices
As far as I know, this is identical to the yDimIndices.
proshade_signed zAxisOrigin
This is the origin position along the z axis.
std::string getRecommendedSymmetryType(void)
This function simply returns the detected recommended symmetry type.
void shiftToRotationCentre(ProSHADE_settings *settings)
Function for shifting map so that its rotation centre is at the centre of the box.
void setPDBMapValues(void)
Function for determining iterator start and stop positions.
void findMapCOM(void)
This function finds the centre of mass of the internal map representation.
void getRealRotFunction(double *rotFunReal, int len)
This function fills the input array with the real rotation function values.
proshade_unsign yDimIndicesOriginal
This is the size of the map cell y dimension in indices.
proshade_unsign getXDim(void)
This function allows access to the map size in indices along the X axis.
std::vector< proshade_double * > * getCyclicAxes(void)
This function allows access to the list of detected cyclic axes.
proshade_double originalPdbTransX
The optimal translation vector as it relates to the original PDB positions (and not the ProSHADE inte...
void rotateMapReciprocalSpace(ProSHADE_settings *settings, proshade_double eulerAlpha, proshade_double eulerBeta, proshade_double eulerGamma)
This function rotates a map based on the given Euler angles.
proshade_signed zTo
This is the final index along the z axis.
void readInStructure(std::string fName, proshade_unsign inputO, ProSHADE_settings *settings, proshade_double *maskArr=nullptr, proshade_unsign maskXDim=0, proshade_unsign maskYDim=0, proshade_unsign maskZDim=0, proshade_double *weightsArr=nullptr, proshade_unsign weigXDim=0, proshade_unsign weigYDim=0, proshade_unsign weigZDim=0)
This function initialises the basic ProSHADE_data variables and reads in a single structure.
proshade_signed xTo
This is the final index along the x axis.
void writeMask(std::string fName, proshade_double *mask)
Function for writing out a mask in MRC MAP format.
proshade_double integrationWeight
The Pearson's c.c. type weighting for the integration.
void mapToSpheres(ProSHADE_settings *settings)
This function converts the internal map onto a set of concentric spheres.
proshade_signed xAxisOrigin
This is the origin position along the x axis.
proshade_double originalPdbTransY
The optimal translation vector as it relates to the original PDB positions (and not the ProSHADE inte...
proshade_double *** rrpMatrices
The energy levels descriptor shell correlation tables.
proshade_signed yTo
This is the final index along the y axis.
proshade_single xDimSize
This is the size of the map cell x dimension in Angstroms.
void computeRotatedSH(void)
This function multiplies the objects spherical harmonics with the Wigner D matrices,...
void prepareFSCFourierMemory(proshade_signed *&cutIndices, fftw_complex *&fCoeffsCut, proshade_signed *noBins, proshade_double **&bindata, proshade_signed *&binCounts, proshade_double *&fscByBin, proshade_single resolution, proshade_signed *cutXDim, proshade_signed *cutYDim, proshade_signed *cutZDim)
This function allocates the memory and makes all preparations required for FSC computation.
void getPredictedTetrahedralSymmetriesList(ProSHADE_settings *settings, std::vector< proshade_double * > *CSymList, proshade_signed *&cutIndices, fftw_complex *&fCoeffsCut, proshade_signed noBins, proshade_double **&bindata, proshade_signed *&binCounts, proshade_double *&fscByBin, proshade_signed xDim, proshade_signed yDim, proshade_signed zDim)
This function predicts a list of all T symmetry axes from the already computed C symmetries list.
proshade_double mapMovFromsChangeY
When the map is translated, the yFrom and yTo values are changed. This variable holds how much they h...
proshade_signed yAxisOrigin
This is the origin position along the y axis.
void computeTranslationMap(ProSHADE_internal_data::ProSHADE_data *obj1)
This function does the computation of the translation map and saves results internally.
proshade_unsign yAxisOrder
This is the order of the y axis.
proshade_double * getRealSphHarmValue(proshade_unsign band, proshade_unsign order, proshade_unsign shell)
This function allows access to the private internal real spherical harmonics values.
proshade_double * internalMap
The internal map data representation, which may be amended as the run progresses.
void writePdb(std::string fName, proshade_double euA=0.0, proshade_double euB=0.0, proshade_double euG=0.0, proshade_double trsX=0.0, proshade_double trsY=0.0, proshade_double trsZ=0.0, proshade_double rotX=0.0, proshade_double rotY=0.0, proshade_double rotZ=0.0, bool firstModel=true)
This function writes out the co-ordinates file with ProSHADE type rotation and translation applied.
proshade_unsign zDimIndicesOriginal
This is the size of the map cell z dimension in indices.
void reportOverlayResults(ProSHADE_settings *settings, std::vector< proshade_double > *rotationCentre, std::vector< proshade_double > *eulerAngles, std::vector< proshade_double > *finalTranslation)
This function reports the results of the overlay mode.
proshade_signed xFrom
This is the starting index along the x axis.
proshade_double computeFSC(ProSHADE_settings *settings, std::vector< proshade_double * > *CSym, size_t symIndex, proshade_signed *&cutIndices, fftw_complex *&fCoeffsCut, proshade_signed noBins, proshade_double **&bindata, proshade_signed *&binCounts, proshade_double *&fscByBin, proshade_signed xDim, proshade_signed yDim, proshade_signed zDim, proshade_unsign rotNumber=1)
This function computes FSC for any given axis in the supplied CSym symmetry axes vector.
void writeMap(std::string fName, std::string title="Created by ProSHADE and written by GEMMI", int mode=2)
Function for writing out the internal structure representation in MRC MAP format.
bool shellBandExists(proshade_unsign shell, proshade_unsign bandVal)
This function checks if particular shell has a particular band.
proshade_unsign getMaxBand(void)
This function returns the maximum band value for the object.
void shiftToBoxCentre(ProSHADE_settings *settings)
Function for shifting map so that its centre of box is at required position.
void translateMap(proshade_double trsX, proshade_double trsY, proshade_double trsZ)
This function simply translates the map by a given number of Angstroms along the three axes.
void invertSHCoefficients(void)
This function computes the shell mapped data from inverting the Spherical Harmonics coefficients.
void writeGemmi(std::string fName, gemmi::Structure gemmiStruct, proshade_double euA=0.0, proshade_double euB=0.0, proshade_double euG=0.0, proshade_double trsX=0.0, proshade_double trsY=0.0, proshade_double trsZ=0.0, proshade_double rotX=0.0, proshade_double rotY=0.0, proshade_double rotZ=0.0, bool firstModel=true)
This function writes out the gemmi::Structure object with ProSHADE type rotation and translation appl...
proshade_unsign getYDim(void)
This function allows access to the map size in indices along the Y axis.
void getImagEMatrixValuesForLM(proshade_signed band, proshade_signed order1, double *eMatsLMImag, int len)
This function fills the input array with the imaginary E matrix values for particular band and order1...
void maskMap(ProSHADE_settings *settings)
Function for computing the map mask using blurring and X IQRs from median.
~ProSHADE_data(void)
Destructor for the ProSHADE_data class.
This class stores all the settings and is passed to the executive classes instead of a multitude of p...
proshade_double xCom
The COM of the map after processing along the X-axis.
void getRealTranslationFunction(double *trsFunReal, int len)
This function fills the input array with the real translation function values.
proshade_double * getImagSphHarmValue(proshade_unsign band, proshade_unsign order, proshade_unsign shell)
This function allows access to the private internal imaginary spherical harmonics values.
proshade_double originalPdbTransZ
The optimal translation vector as it relates to the original PDB positions (and not the ProSHADE inte...
void writeOutOverlayFiles(ProSHADE_settings *settings, proshade_double eulA, proshade_double eulB, proshade_double eulG, std::vector< proshade_double > *rotCentre, std::vector< proshade_double > *ultimateTranslation)
This function writes out the rotated map, co-ordinates and transformation JSON file.
void getImagRotFunction(double *rotFunImag, int len)
This function fills the input array with the imaginary rotation function values.
void processInternalMap(ProSHADE_settings *settings)
This function simply clusters several other functions which should be called together.
proshade_complex * getTranslationFnPointer(void)
This function allows access to the translation function through a pointer.
proshade_double zCom
The COM of the map after processing along the Z-axis.
std::vector< proshade_double > rotateMapRealSpaceInPlace(proshade_double eulA, proshade_double eulB, proshade_double eulG)
This function rotates a map based on the given Euler angles in place.
std::vector< std::vector< proshade_double * > > dihedralSymmetries
This is where the detected dihedral ("D") symmetries will be kept.
std::vector< proshade_single > spherePos
Vector of sphere radii from the centre of the map.
proshade_unsign getZDim(void)
This function allows access to the map size in indices along the Z axis.
proshade_single getSpherePosValue(proshade_unsign shell)
This function allows access to sphere positions.
void convertRotationFunction(ProSHADE_settings *settings)
This function converts the self-rotation function of this structure to angle-axis representation.
proshade_double getRRPValue(proshade_unsign band, proshade_unsign sh1, proshade_unsign sh2)
This function allows access to the priva internal RRP matrices.
proshade_double mapMovFromsChangeZ
When the map is translated, the zFrom and zTo values are changed. This variable holds how much they h...
proshade_unsign recommendedSymmetryFold
The fold of the recommended symmetry C or D type, 0 otherwise.
std::vector< proshade_double * > recommendedSymmetryValues
The axes and other info of the recommended symmetry for the structure.
proshade_double getMapValue(proshade_unsign pos)
This function returns the internal map representation value of a particular array position.
void getImagSO3Coeffs(double *so3CoefsImag, int len)
This function fills the input array with the imaginary SO(3) coefficient values.
proshade_signed * getYToPtr(void)
This function allows access to the map last position along the Y axis.
proshade_complex *** wignerMatrices
These matrices are computed for a particular rotation to be done in spherical harmonics.
ProSHADE_internal_io::InputType fileType
This is the type of the input file.
std::vector< proshade_double > rotateMapRealSpace(proshade_double axX, proshade_double axY, proshade_double axZ, proshade_double axAng, proshade_double *&map)
This function rotates a map based on the given angle-axis rotation.
proshade_signed * getZAxisOrigin(void)
This function allows access to the map Z axis origin value.
proshade_single yDimSizeOriginal
This is the size of the map cell y dimension in Angstroms.
proshade_unsign zAxisOrder
This is the order of the z axis.
void reSampleMap(ProSHADE_settings *settings)
This function changes the internal map sampling to conform to particular resolution value.
void rotateFourierCoeffs(proshade_double axX, proshade_double axY, proshade_double axZ, proshade_double axAng, fftw_complex *&coeffs, fftw_complex *&rotCoeffs, proshade_signed xDim, proshade_signed yDim, proshade_signed zDim)
This function takes the rotation in terms of angle-axis representation and returns the input Fourier ...
proshade_unsign maxEMatDim
The band (l) value for E matrix (i.e. the smallest of the two bands).
proshade_signed * getXToPtr(void)
This function allows access to the map last position along the X axis.
proshade_signed * getYFromPtr(void)
This function allows access to the map start along the Y axis.
proshade_double mapCOMProcessChangeY
The change in Y axis between the creation of the structure (originalMapYCom) and just before rotation...
proshade_single bAngle
This is the angle b of the map cell in degrees.
proshade_signed yAxisOriginOriginal
This is the origin position along the y axis.
proshade_single aAngle
This is the angle a of the map cell in degrees.
This namespace contains the structure and functions required for data reading and storing their deriv...
void reportSymmetryResults(ProSHADE_settings *settings)
This function prints the report for symmetry detection.
void figureIndexStartStop(void)
Function for determining iterator start and stop positions.
bool isEmpty
This variable stated whether the class contains any information.
void reportCurrentSymmetryResults(ProSHADE_settings *settings, proshade_double threshold, proshade_signed *&cutIndices, fftw_complex *&fCoeffsCut, proshade_signed noBins, proshade_double **&bindata, proshade_signed *&binCounts, proshade_double *&fscByBin, proshade_signed cutXDim, proshade_signed cutYDim, proshade_signed cutZDim)
This function prints a single line of the symmetry detection results for a particular threshold.
void setEMatrixValue(int band, int order1, int order2, proshade_complex val)
This function allows setting the E matrix value.
proshade_complex * getInvSO3Coeffs(void)
This function allows access to the inverse SO(3) coefficients array.
proshade_unsign getEMatDim(void)
This function allows access to the maximum band for the E matrix.
This header file declares functions required for 3D map peak searching.
std::string recommendedSymmetryType
The symmetry type that ProSHADE finds the best fitting for the structure. Possible values are "" for ...
proshade_signed yFrom
This is the starting index along the y axis.
void allocateEMatrices(proshade_unsign band, proshade_single oversamplingRatio)
This function allocates the required memory for the E matrices.
proshade_signed * getYAxisOrigin(void)
This function allows access to the map Y axis origin value.
std::vector< proshade_double * > getCyclicSymmetriesListFromAngleAxis(ProSHADE_settings *settings)
This function obtains a list of all C symmetries from the angle-axis space mapped rotation function v...
proshade_double originalPdbRotCenZ
The centre of rotation as it relates to the original PDB positions (and not the ProSHADE internal map...
ProSHADE_data()
Constructor for getting empty ProSHADE_data class.
proshade_unsign getNoRecommendedSymmetryAxes(void)
This function returns the number of detected recommended symmetry axes.
void reportSymmetryResultsList(ProSHADE_settings *settings)
This function takes prints the report for symmetry detection using multiple thresholds....
proshade_unsign zDimIndices
This is the size of the map cell z dimension in indices.
std::vector< std::vector< proshade_double > > computeGroupElementsForGroup(proshade_double xAx, proshade_double yAx, proshade_double zAx, proshade_signed fold)
This function computes the group elements as rotation matrices (except for the identity element) for ...
proshade_unsign zGridIndices
As far as I know, this is identical to the zDimIndices.
std::vector< proshade_double * > tetrahedralSymmetries
This is where the detected tetrahedral ("T") symmetries will be kept.
void getDihedralSymmetriesList(ProSHADE_settings *settings, std::vector< proshade_double * > *CSymList)
This function obtains a list of all D symmetries from already computed C symmetries list.
ProSHADE_internal_spheres::ProSHADE_sphere ** spheres
The set of concentric spheres to which the intermal density map has been projected.
void normaliseMap(ProSHADE_settings *settings)
Function for normalising the map values to mean 0 and sd 1..
std::vector< proshade_double * > findRequestedCSymmetryFromAngleAxis(ProSHADE_settings *settings, proshade_unsign fold, proshade_double *peakThres)
This function searches the angle-axis representation of the rotation function for a cyclic point grou...
void detectSymmetryFromAngleAxisSpace(ProSHADE_settings *settings)
This function runs the symmetry detection algorithms on this structure using the angle-axis space and...
void allocateSO3CoeffsSpace(proshade_unsign band)
This function allocates the memory for the SO(3) coefficients and the inverse for that calling object...
proshade_double originalPdbRotCenX
The centre of rotation as it relates to the original PDB positions (and not the ProSHADE internal map...
proshade_double originalMapZCom
The COM of the first map to be loaded/computed without any furhter changes being reflacted along the ...
void getRealSO3Coeffs(double *so3CoefsReal, int len)
This function fills the input array with the real SO(3) coefficient values.
proshade_single yDimSize
This is the size of the map cell y dimension in Angstroms.
void removePhaseInormation(ProSHADE_settings *settings)
This function removes phase from the map, effectively converting it to Patterson map.
void getPredictedOctahedralSymmetriesList(ProSHADE_settings *settings, std::vector< proshade_double * > *CSymList, proshade_signed *&cutIndices, fftw_complex *&fCoeffsCut, proshade_signed noBins, proshade_double **&bindata, proshade_signed *&binCounts, proshade_double *&fscByBin, proshade_signed xDim, proshade_signed yDim, proshade_signed zDim)
This function predicts a list of all O symmetry axes from the already computed C symmetries list.
void zeroPaddToDims(proshade_unsign xDim, proshade_unsign yDim, proshade_unsign zDim)
This function changes the size of a structure to fit the supplied new limits.
void addExtraSpace(ProSHADE_settings *settings)
This function increases the size of the map so that it can add empty space around it.
std::vector< proshade_double * > octahedralSymmetries
This is where the detected octahedral ("O") symmetries will be kept.
void getRealEMatrixValuesForLM(proshade_signed band, proshade_signed order1, double *eMatsLMReal, int len)
This function fills the input array with the real E matrix values for particular band and order1 (l a...
void readInPDB(ProSHADE_settings *settings)
Function for reading pdb data.
proshade_signed * getZFromPtr(void)
This function allows access to the map start along the Z axis.
proshade_signed * getZToPtr(void)
This function allows access to the map last position along the Z axis.
int so3CoeffsArrayIndex(proshade_signed order1, proshade_signed order2, proshade_signed band)
This function gets the SO(3) coefficients array index for a particular so(3) band,...