![]() |
ProSHADE
0.7.6.1 (AUG 2021)
Protein Shape Detection
|
This namespace contains the internal functions for manipulating maps already present in the internal structures. More...
Functions | |
proshade_signed | myRound (proshade_double x) |
Calls the appropriate version of round function depending on compiler version. More... | |
proshade_signed | myRound (proshade_single x) |
Calls the appropriate version of round function depending on compiler version. More... | |
void | determinePDBRanges (gemmi::Structure pdbFile, proshade_single *xFrom, proshade_single *xTo, proshade_single *yFrom, proshade_single *yTo, proshade_single *zFrom, proshade_single *zTo, bool firstModel) |
Function for finding the PDB file ranges. More... | |
void | findPDBCOMValues (gemmi::Structure pdbFile, proshade_double *xCom, proshade_double *yCom, proshade_double *zCom, bool firstModel) |
This function finds the Centre of Mass for the co-ordinate file. More... | |
void | findMAPCOMValues (proshade_double *map, proshade_double *xCom, proshade_double *yCom, proshade_double *zCom, proshade_single xAngs, proshade_single yAngs, proshade_single zAngs, proshade_signed xFrom, proshade_signed xTo, proshade_signed yFrom, proshade_signed yTo, proshade_signed zFrom, proshade_signed zTo) |
This function finds the Centre of Mass for a map. More... | |
void | rotatePDBCoordinates (gemmi::Structure *pdbFile, proshade_double euA, proshade_double euB, proshade_double euG, proshade_double xCom, proshade_double yCom, proshade_double zCom, bool firstModel) |
Function for rotating the PDB file co-ordinates by Euler angles. More... | |
void | translatePDBCoordinates (gemmi::Structure *pdbFile, proshade_double transX, proshade_double transY, proshade_double transZ, bool firstModel) |
Function for translating the PDB file co-ordinates by given distances in Angstroms. More... | |
void | changePDBBFactors (gemmi::Structure *pdbFile, proshade_double newBFactorValue, bool firstModel) |
Function for changing the PDB B-factor values to a specific single value. More... | |
void | removeWaters (gemmi::Structure *pdbFile, bool firstModel) |
This function removed all waters from PDB input file. More... | |
void | movePDBForMapCalc (gemmi::Structure *pdbFile, proshade_single xMov, proshade_single yMov, proshade_single zMov, bool firstModel) |
Function for moving co-ordinate atoms to better suit theoretical map computation. More... | |
void | generateMapFromPDB (gemmi::Structure pdbFile, proshade_double *&map, proshade_single requestedResolution, proshade_single xCell, proshade_single yCell, proshade_single zCell, proshade_signed *xTo, proshade_signed *yTo, proshade_signed *zTo, bool forceP1, bool firstModel) |
This function generates a theoretical map from co-ordinate input files. More... | |
void | moveMapByIndices (proshade_single *xMov, proshade_single *yMov, proshade_single *zMov, proshade_single xAngs, proshade_single yAngs, proshade_single zAngs, proshade_signed *xFrom, proshade_signed *xTo, proshade_signed *yFrom, proshade_signed *yTo, proshade_signed *zFrom, proshade_signed *zTo, proshade_signed *xOrigin, proshade_signed *yOrigin, proshade_signed *zOrigin) |
Function for moving map back to original PDB location by changing the indices. More... | |
void | moveMapByFourier (proshade_double *&map, proshade_single xMov, proshade_single yMov, proshade_single zMov, proshade_single xAngs, proshade_single yAngs, proshade_single zAngs, proshade_signed xDim, proshade_signed yDim, proshade_signed zDim) |
Function for moving map back to original PDB location by using Fourier transformation. More... | |
void | blurSharpenMap (proshade_double *&map, proshade_double *&maskedMap, proshade_unsign xDimS, proshade_unsign yDimS, proshade_unsign zDimS, proshade_single xAngs, proshade_single yAngs, proshade_single zAngs, proshade_single blurringFactor) |
Function for blurring/sharpening maps. More... | |
void | getMaskFromBlurr (proshade_double *&blurMap, proshade_double *&outMap, proshade_unsign xDimS, proshade_unsign yDimS, proshade_unsign zDimS, proshade_single noIQRs) |
Function for computing mask from blurred map. More... | |
void | getNonZeroBounds (proshade_double *map, proshade_signed xDim, proshade_signed yDim, proshade_signed zDim, proshade_signed *&ret) |
Function for finding the map boundaries enclosing positive only values. More... | |
void | addExtraBoundSpace (proshade_unsign xDim, proshade_unsign yDim, proshade_unsign zDim, proshade_single xAngs, proshade_single yAngs, proshade_single zAngs, proshade_signed *&bounds, proshade_single extraSpace) |
This function takes a set of bounds and adds a given number of Angstroms to them. More... | |
void | reSampleMapToResolutionTrilinear (proshade_double *&map, proshade_single resolution, proshade_unsign xDimS, proshade_unsign yDimS, proshade_unsign zDimS, proshade_single xAngs, proshade_single yAngs, proshade_single zAngs, proshade_single *&corrs) |
This function re-samples a map to conform to given resolution using tri-linear interpolation. More... | |
void | reSampleMapToResolutionFourier (proshade_double *&map, proshade_single resolution, proshade_unsign xDimS, proshade_unsign yDimS, proshade_unsign zDimS, proshade_single xAngs, proshade_single yAngs, proshade_single zAngs, proshade_single *&corrs) |
This function re-samples a map to conform to given resolution using Fourier. More... | |
void | allocateResolutionFourierMemory (fftw_complex *&origMap, fftw_complex *&fCoeffs, fftw_complex *&newFCoeffs, fftw_complex *&newMap, fftw_plan &planForwardFourier, fftw_plan &planBackwardRescaledFourier, proshade_unsign xDimOld, proshade_unsign yDimOld, proshade_unsign zDimOld, proshade_unsign xDimNew, proshade_unsign yDimNew, proshade_unsign zDimNew) |
This function allocates and checks the allocatio of the memory required by the Fourier resampling. More... | |
void | releaseResolutionFourierMemory (fftw_complex *&origMap, fftw_complex *&fCoeffs, fftw_complex *&newFCoeffs, fftw_complex *&newMap, fftw_plan &planForwardFourier, fftw_plan &planBackwardRescaledFourier) |
This function releases the memory required by the Fourier resampling. More... | |
void | changeFourierOrder (fftw_complex *&fCoeffs, proshade_signed xDim, proshade_signed yDim, proshade_signed zDim, bool negativeFirst) |
This function changes the order of Fourier coefficients in a 3D array between positive first (default) and negative first (mass centered at xMax/2, yMax/2, zMax/2 instead of 0,0,0) More... | |
void | removeMapPhase (fftw_complex *&mapCoeffs, proshade_unsign xDim, proshade_unsign yDim, proshade_unsign zDim) |
This function removes the phase from reciprocal (frequency) map. More... | |
void | getFakeHalfMap (proshade_double *&map, proshade_double *&fakeHalfMap, proshade_unsign xDimS, proshade_unsign yDimS, proshade_unsign zDimS, proshade_signed fakeMapKernel) |
Function for creating "fake" half-maps. More... | |
void | getCorrelationMapMask (proshade_double *&map, proshade_double *&fakeHalfMap, proshade_double *&correlationMask, proshade_unsign xDimS, proshade_unsign yDimS, proshade_unsign zDimS, proshade_signed corrMaskKernel) |
Function for creating the correlation mask. More... | |
proshade_single | getIndicesFromAngstroms (proshade_unsign xDim, proshade_unsign yDim, proshade_unsign zDim, proshade_single xAngs, proshade_single yAngs, proshade_single zAngs, proshade_single dist) |
This function converts distance in Angstroms to distance in map indices. More... | |
void | connectMaskBlobs (proshade_double *&mask, proshade_signed xDim, proshade_signed yDim, proshade_signed zDim, proshade_single xAngs, proshade_single yAngs, proshade_single zAngs, proshade_single maskThres) |
This function connects blobs in mask. More... | |
void | beautifyBoundaries (proshade_signed *&bounds, proshade_unsign xDim, proshade_unsign yDim, proshade_unsign zDim, proshade_signed boundsDiffThres) |
Function for modifying boundaries to a mathematically more pleasant values. More... | |
proshade_signed | betterClosePrimeFactors (proshade_signed fromRange, proshade_signed toRange) |
Function for finding close numbers with better prime factors. More... | |
void | distributeSpaceToBoundaries (proshade_signed &minBound, proshade_signed &maxBound, proshade_signed oldBoundRange, proshade_signed newBoundRange) |
Function for adding space to boundaries within the map confines. More... | |
void | copyMapByBounds (proshade_signed xFrom, proshade_signed xTo, proshade_signed yFrom, proshade_signed yTo, proshade_signed zFrom, proshade_signed zTo, proshade_signed origXFrom, proshade_signed origYFrom, proshade_signed origZFrom, proshade_unsign yDimIndices, proshade_unsign zDimIndices, proshade_unsign origXDimIndices, proshade_unsign origYDimIndices, proshade_unsign origZDimIndices, proshade_double *&newMap, proshade_double *origMap) |
This function copies an old map to a new map with different boundaries. More... | |
This namespace contains the internal functions for manipulating maps already present in the internal structures.
The ProSHADE_internal_mapManip namespace contains helper functions for map manipulation and processing. However, these functions do make minimum assumptions about the map internal organisation and variables and simply perform tasks on generic variables. None of these functions should be used directly be the user.
void ProSHADE_internal_mapManip::addExtraBoundSpace | ( | proshade_unsign | xDim, |
proshade_unsign | yDim, | ||
proshade_unsign | zDim, | ||
proshade_single | xAngs, | ||
proshade_single | yAngs, | ||
proshade_single | zAngs, | ||
proshade_signed *& | bounds, | ||
proshade_single | extraSpace | ||
) |
This function takes a set of bounds and adds a given number of Angstroms to them.
This function adds an extra distance to a set of bounds as given by the number of angstroms to be added. If the resulting bounds are larger than the maximum map bounds, the maxima and minima are returned instead, thus partially ignoring the requested extra space. !!! If the resulting index ranges would be odd, extra one idice is added to make them even. !!!
[in] | xDim | The number of indices along the x axis of the map. |
[in] | yDim | The number of indices along the y axis of the map. |
[in] | zDim | The number of indices along the z axis of the map. |
[in] | xAngs | The size of the x dimension of the map in angstroms. |
[in] | yAngs | The size of the y dimension of the map in angstroms. |
[in] | zAngs | The size of the z dimension of the map in angstroms. |
[in] | bounds | A proshade_unsign pointer reference to array of 6 which has the current bounds (0 = minX; 1 = maxX; 2 = minY; 3 = maxY; 4 - minZ; 5 = maxZ). |
[in] | extraSpace | The number of angstroms to be added to each dimension (both to start and to end). |
Definition at line 1139 of file ProSHADE_mapManip.cpp.
void ProSHADE_internal_mapManip::allocateResolutionFourierMemory | ( | fftw_complex *& | origMap, |
fftw_complex *& | fCoeffs, | ||
fftw_complex *& | newFCoeffs, | ||
fftw_complex *& | newMap, | ||
fftw_plan & | planForwardFourier, | ||
fftw_plan & | planBackwardRescaledFourier, | ||
proshade_unsign | xDimOld, | ||
proshade_unsign | yDimOld, | ||
proshade_unsign | zDimOld, | ||
proshade_unsign | xDimNew, | ||
proshade_unsign | yDimNew, | ||
proshade_unsign | zDimNew | ||
) |
This function allocates and checks the allocatio of the memory required by the Fourier resampling.
This function allocates the memory required for the Fourier space re-sampling of density maps. It allocates the original map and original map coefficients arrays (these need to be of the fftw_complex type) as well as the re-sampled map and re-sampled map coefficients. It then also proceeds to check the memory allocation and creating the FFTW transform plans for both, the forward and re-sampled backward operations.
[in] | origMap | A Reference pointer to an array where the original map data will be stored. |
[in] | fCoeffs | A Reference pointer to an array where the original map Fourier coefficients data will be stored. |
[in] | newFCoeffs | A Reference pointer to an array where the re-sampled map Fourier coefficients data will be stored. |
[in] | newMap | A Reference pointer to an array where the re-sampled map data will be stored. |
[in] | planForwardFourier | FFTW_plan which will compute the original map to Fourier coefficients transform. |
[in] | planBackwardRescaledFourier | FFTW_plan which will compute the re-sampled Fourier coefficients to re-sampled map transform. |
[in] | xDimOld | The number of indices along the x-axis of the original map. |
[in] | yDimOld | The number of indices along the y-axis of the original map. |
[in] | zDimOld | The number of indices along the z-axis of the original map. |
[in] | xDimNew | The number of indices along the x-axis of the re-sampled map. |
[in] | yDimNew | The number of indices along the y-axis of the re-sampled map. |
[in] | zDimNew | The number of indices along the z-axis of the re-sampled map. |
Definition at line 1496 of file ProSHADE_mapManip.cpp.
void ProSHADE_internal_mapManip::beautifyBoundaries | ( | proshade_signed *& | bounds, |
proshade_unsign | xDim, | ||
proshade_unsign | yDim, | ||
proshade_unsign | zDim, | ||
proshade_signed | boundsDiffThres | ||
) |
Function for modifying boundaries to a mathematically more pleasant values.
This function does two things. Firstly, it attempts to find even boundaries which have a small sum of the prime factors - i.e. this is good for further processing. Secondly, it attempts to make boundaries with sizes within some threshold the same. This is purely human thing :-).
[in] | bounds | A proshade_signed pointer to array of 6 in which the current bounds are (0 = minX; 1 = maxX; 2 = minY; 3 = maxY; 4 - minZ; 5 = maxZ). |
[in] | xDim | The number of indices along the x axis of the map. |
[in] | yDim | The number of indices along the y axis of the map. |
[in] | zDim | The number of indices along the z axis of the map. |
[in] | boundsDiffThres | Number of indices by which boudaries in different dimensions can differ and still be made the same. |
Definition at line 1916 of file ProSHADE_mapManip.cpp.
proshade_signed ProSHADE_internal_mapManip::betterClosePrimeFactors | ( | proshade_signed | fromRange, |
proshade_signed | toRange | ||
) |
Function for finding close numbers with better prime factors.
This function takes a range of integer numbers and proceeds to find if a close number in the range to the range start does not have better (smaller sum of) prime factors, taking into account the distance of this better number to the starting position. This is useful for finding boundaries with nice prime number sizes.
[in] | fromRange | Which number the search will be from and possibilities compared to. |
[in] | toRange | The maximum number to which the posibilities should be searched. |
[out] | X | The value in the range which has better prime factors, or the first value if none other has. |
Definition at line 1986 of file ProSHADE_mapManip.cpp.
void ProSHADE_internal_mapManip::blurSharpenMap | ( | proshade_double *& | map, |
proshade_double *& | blurredMap, | ||
proshade_unsign | xDimS, | ||
proshade_unsign | yDimS, | ||
proshade_unsign | zDimS, | ||
proshade_single | xAngs, | ||
proshade_single | yAngs, | ||
proshade_single | zAngs, | ||
proshade_single | blurringFactor | ||
) |
Function for blurring/sharpening maps.
This function takes the internal map representation information from the parameters as well as the internal map itself and proceeds to compute blue/sharpen the map by changing its overall B-factors by adding a specific value. If this value is negative, then the map is sharpened, while if the value is positive, the map is blurred.
[in] | map | A Reference Pointer to the map which should be blurred/sharpened. |
[in] | blurredMap | A Reference Pointer to the variable which will store the modified map. |
[in] | xDimS | The number of indices along the x axis of the map. |
[in] | yDimS | The number of indices along the y axis of the map. |
[in] | zDimS | The number of indices along the z axis of the map. |
[in] | xAngs | The size of the x dimension of the map in angstroms. |
[in] | yAngs | The size of the y dimension of the map in angstroms. |
[in] | zAngs | The size of the z dimension of the map in angstroms. |
[in] | blurringFactor | The amount of B-factor change that should be applied to the map. (I.e. increasing the map overall B-factors by 20 will be done by supplying 20 as a value for this parameter). |
Definition at line 932 of file ProSHADE_mapManip.cpp.
void ProSHADE_internal_mapManip::changeFourierOrder | ( | fftw_complex *& | fCoeffs, |
proshade_signed | xDim, | ||
proshade_signed | yDim, | ||
proshade_signed | zDim, | ||
bool | negativeFirst | ||
) |
This function changes the order of Fourier coefficients in a 3D array between positive first (default) and negative first (mass centered at xMax/2, yMax/2, zMax/2 instead of 0,0,0)
This function firstly determines the start and end of the positive and negative Fourier coefficients for all three axes (it assumes 3D array); this works for both odd and even axis sizes. To do this properly, in the case of odd axis sizes, the function needs to know whether we are changing the order from FFTW defaul positive first, or if we are changing the other way around - the last parameter serves the purpose of passing this information. Finally, the function then switches the order of the Fourier coefficients as requested using a temporary array that it allocates and deletes within itself.
[in] | fCoeffs | A Reference pointer to an array where the Fourier coefficients for re-ordering are stored. |
[in] | xDim | The size of the x-axis dimension of the fCoeffs array. |
[in] | yDim | The size of the x-axis dimension of the fCoeffs array. |
[in] | zDim | The size of the x-axis dimension of the fCoeffs array. |
[in] | negativeFirst | Should the coefficients be stored negarive first (TRUE), or are we reversing already re-ordered array (FALSE)? |
Definition at line 1559 of file ProSHADE_mapManip.cpp.
void ProSHADE_internal_mapManip::changePDBBFactors | ( | gemmi::Structure * | pdbFile, |
proshade_double | newBFactorValue, | ||
bool | firstModel | ||
) |
Function for changing the PDB B-factor values to a specific single value.
This function does a quick read-through the PDB file and changes all the atom B-factor values to a single, pre-set value. This is important for PDB files which have all atoms with B-factor value 0, as these then produce bad maps, which fail further processing by ProSHADE.
[in] | pdbFile | A pointer to gemmi::Structure object read in from the input file. |
[in] | newBFactorValue | The value with which all atom B-factor values should be replaced. |
[in] | firstModel | Should only the first, or all models be used? |
Definition at line 444 of file ProSHADE_mapManip.cpp.
void ProSHADE_internal_mapManip::connectMaskBlobs | ( | proshade_double *& | mask, |
proshade_signed | xDim, | ||
proshade_signed | yDim, | ||
proshade_signed | zDim, | ||
proshade_single | xAngs, | ||
proshade_single | yAngs, | ||
proshade_single | zAngs, | ||
proshade_single | maskThres | ||
) |
This function connects blobs in mask.
[in] | mask | A pointer reference to mask map in which blobs should be connected. |
[in] | xDim | The number of map indices along the x-axis. |
[in] | yDim | The number of map indices along the y-axis. |
[in] | zDim | The number of map indices along the z-axis. |
[in] | xAngs | The map size in Angstroms along the x-axis. |
[in] | yAngs | The map size in Angstroms along the y-axis. |
[in] | zAngs | The map size in Angstroms along the z-axis. |
[in] | maskThres | The threshold which will be used for applying mask. |
Definition at line 1839 of file ProSHADE_mapManip.cpp.
void ProSHADE_internal_mapManip::copyMapByBounds | ( | proshade_signed | xFrom, |
proshade_signed | xTo, | ||
proshade_signed | yFrom, | ||
proshade_signed | yTo, | ||
proshade_signed | zFrom, | ||
proshade_signed | zTo, | ||
proshade_signed | origXFrom, | ||
proshade_signed | origYFrom, | ||
proshade_signed | origZFrom, | ||
proshade_unsign | yDimIndices, | ||
proshade_unsign | zDimIndices, | ||
proshade_unsign | origXDimIndices, | ||
proshade_unsign | origYDimIndices, | ||
proshade_unsign | origZDimIndices, | ||
proshade_double *& | newMap, | ||
proshade_double * | origMap | ||
) |
This function copies an old map to a new map with different boundaries.
This function takes old and new structure details and the original structure map. It then proceed to iterate through the new structure map, checking if this particular point is inside the original structure's map and if it is, it will copy the value of this point from the original to the new map. Otherwise, it will set the new map's point value to 0.0. This is used when creating a new structure from an old structure given new bounds, which can be larger or smaller than the originals.
[in] | xFrom | The starting x-axis index of the new map. |
[in] | xTo | The last x-axis index of the new map. |
[in] | yFrom | The starting y-axis index of the new map. |
[in] | yTo | The last y-axis index of the new map. |
[in] | zFrom | The starting z-axis index of the new map. |
[in] | zTo | The last z-axis index of the new map. |
[in] | origXFrom | The starting x-axis index of the original (copied) map. |
[in] | origYFrom | The starting y-axis index of the original (copied) map. |
[in] | origZFrom | The starting z-axis index of the original (copied) map. |
[in] | yDimIndices | The size of the y-axis dimension in the new map. |
[in] | zDimIndices | The size of the z-axis dimension in the new map. |
[in] | origXDimIndices | The size of the x-axis dimension in the old map. |
[in] | origYDimIndices | The size of the y-axis dimension in the old map. |
[in] | origZDimIndices | The size of the z-axis dimension in the old map. |
[in] | newMap | Pointer reference to the array where new map will be saved to. |
[in] | origMap | Pointer to the array where the original map can be accessed. |
Definition at line 2082 of file ProSHADE_mapManip.cpp.
void ProSHADE_internal_mapManip::determinePDBRanges | ( | gemmi::Structure | pdbFile, |
proshade_single * | xFrom, | ||
proshade_single * | xTo, | ||
proshade_single * | yFrom, | ||
proshade_single * | yTo, | ||
proshade_single * | zFrom, | ||
proshade_single * | zTo, | ||
bool | firstModel | ||
) |
Function for finding the PDB file ranges.
This function does a quick read-through the PDB file and reports the x, y and z to and from values. This is used to determine if these need to be changed for proper theoretical map computation operation.
[in] | pdbFile | A gemmi::Structure object read in from the input file. |
[in] | xFrom | Address to a variable to save the x axis minimal atom position. |
[in] | xTo | Address to a variable to save the x axis maximum atom position. |
[in] | yFrom | Address to a variable to save the y axis minimal atom position. |
[in] | yTo | Address to a variable to save the y axis maximum atom position. |
[in] | zFrom | Address to a variable to save the z axis minimal atom position. |
[in] | zTo | Address to a variable to save the z axis maximum atom position. |
[in] | firstModel | Should only the first, or all models be used? |
Definition at line 72 of file ProSHADE_mapManip.cpp.
void ProSHADE_internal_mapManip::distributeSpaceToBoundaries | ( | proshade_signed & | minBound, |
proshade_signed & | maxBound, | ||
proshade_signed | oldBoundRange, | ||
proshade_signed | newBoundRange | ||
) |
Function for adding space to boundaries within the map confines.
This function takes the current boudaries and attempts to add the extra space specified by the parameters 3 and 4 equally to both the top and bottom boundaries, making sure the bottom boundary does not go belowe 0 and the top boundary does not exceed the final parameter.
[in] | minBound | Reference to where the lower boundary value is held. |
[in] | maxBound | Reference to where the upper boundary value is held. |
[in] | oldBoundRange | The range between the original boundaries. |
[in] | newBoundRange | The range which should be between the new boundaries. |
Definition at line 2032 of file ProSHADE_mapManip.cpp.
void ProSHADE_internal_mapManip::findMAPCOMValues | ( | proshade_double * | map, |
proshade_double * | xCom, | ||
proshade_double * | yCom, | ||
proshade_double * | zCom, | ||
proshade_single | xAngs, | ||
proshade_single | yAngs, | ||
proshade_single | zAngs, | ||
proshade_signed | xFrom, | ||
proshade_signed | xTo, | ||
proshade_signed | yFrom, | ||
proshade_signed | yTo, | ||
proshade_signed | zFrom, | ||
proshade_signed | zTo | ||
) |
This function finds the Centre of Mass for a map.
This function takes the proshade map object and procceds to compute the Centre of Mass of this object in Angstroms in real space.
[in] | map | A pointer to prshade density map. |
[in] | xCom | A pointer to proshade_double variable where the MAP file COM along the X-axis will be saved. |
[in] | yCom | A pointer to proshade_double variable where the MAP file COM along the Y-axis will be saved. |
[in] | zCom | A pointer to proshade_double variable where the MAP file COM along the Z-axis will be saved. |
[in] | xAngs | The map size in Angstroms along the X axis. |
[in] | yAngs | The map size in Angstroms along the Y axis. |
[in] | zAngs | The map size in Angstroms along the Z axis. |
[in] | xFrom | The initial index of the x dimension of the map. |
[in] | xTo | The terminal index of the x dimension of the map. |
[in] | yFrom | The initial index of the y dimension of the map. |
[in] | yTo | The terminal index of the y dimension of the map. |
[in] | zFrom | The initial index of the z dimension of the map. |
[in] | zTo | The terminal index of the z dimension of the map. |
Definition at line 240 of file ProSHADE_mapManip.cpp.
void ProSHADE_internal_mapManip::findPDBCOMValues | ( | gemmi::Structure | pdbFile, |
proshade_double * | xCom, | ||
proshade_double * | yCom, | ||
proshade_double * | zCom, | ||
bool | firstModel | ||
) |
This function finds the Centre of Mass for the co-ordinate file.
This function takes the gemmi::Structure object read from a co-odinate file and procceds to compute the Centre of Mass of this object.
[in] | pdbFile | A gemmi::Structure object read in from input file. |
[in] | xCom | A pointer to proshade_double variable where the PDB file COM along the X-axis will be saved. |
[in] | yCom | A pointer to proshade_double variable where the PDB file COM along the Y-axis will be saved. |
[in] | zCom | A pointer to proshade_double variable where the PDB file COM along the Z-axis will be saved. |
[in] | firstModel | Should only the first, or all models be used? |
Definition at line 157 of file ProSHADE_mapManip.cpp.
void ProSHADE_internal_mapManip::generateMapFromPDB | ( | gemmi::Structure | pdbFile, |
proshade_double *& | map, | ||
proshade_single | requestedResolution, | ||
proshade_single | xCell, | ||
proshade_single | yCell, | ||
proshade_single | zCell, | ||
proshade_signed * | xTo, | ||
proshade_signed * | yTo, | ||
proshade_signed * | zTo, | ||
bool | forceP1, | ||
bool | firstModel | ||
) |
This function generates a theoretical map from co-ordinate input files.
This function makes use of the Gemmi internal functionality to firstly create a properly sized cell object, then to check the input co-ordinate file for containing known elements as well as elements for which Gemmi knows the form factors. Then, this function proceeds to compute the F's using Cromer & Libermann method (from Gemmi) to finally compute the theoretical map using these. This map is then copied to the variable supplied in the second argument, presumably the ProSHADE internal map variable.
[in] | pdbFile | A gemmi::Structure object read in from the input file. |
[in] | map | Pointer reference to a variable to save the map data. |
[in] | requestedResolution | Map resolution to which the map should be computed. |
[in] | xCell | The size of the map cell to be created. |
[in] | yCell | The size of the map cell to be created. |
[in] | zCell | The size of the map cell to be created. |
[in] | xTo | Pointer to variable where the map size along the x-axis in indices will be saved. |
[in] | yTo | Pointer to variable where the map size along the y-axis in indices will be saved. |
[in] | zTo | Pointer to variable where the map size along the z-axis in indices will be saved. |
[in] | forceP1 | Should the P1 spacegroup be forced? |
[in] | firstModel | Should only the first, or all models be used? |
Definition at line 644 of file ProSHADE_mapManip.cpp.
void ProSHADE_internal_mapManip::getCorrelationMapMask | ( | proshade_double *& | map, |
proshade_double *& | fakeHalfMap, | ||
proshade_double *& | correlationMask, | ||
proshade_unsign | xDimS, | ||
proshade_unsign | yDimS, | ||
proshade_unsign | zDimS, | ||
proshade_signed | corrMaskKernel | ||
) |
Function for creating the correlation mask.
This function is currently not used and should probably be deleted. The proper implementation of this masking approach should be available in the EMDA software.
[in] | map | A Reference Pointer to the map which should be blurred/sharpened. |
[in] | blurredMap | A Reference Pointer to the variable which stores the fake half-map. |
[in] | correlationMask | A Reference Pointer to empty map where the mask will be saved. |
[in] | xDimS | The number of indices along the x axis of the map. |
[in] | yDimS | The number of indices along the y axis of the map. |
[in] | zDimS | The number of indices along the z axis of the map. |
[in] | corrMaskKernel | The amount of neighbours in any direction whose correlation is to be used to get the current points correlation. |
Definition at line 1745 of file ProSHADE_mapManip.cpp.
void ProSHADE_internal_mapManip::getFakeHalfMap | ( | proshade_double *& | map, |
proshade_double *& | fakeHalfMap, | ||
proshade_unsign | xDimS, | ||
proshade_unsign | yDimS, | ||
proshade_unsign | zDimS, | ||
proshade_signed | fakeMapKernel | ||
) |
Function for creating "fake" half-maps.
This function takes the internal map and an empty map array and proceeds to find all neighbours within the kernel distance to each map points. It then computes the average of these neighbours and saves this average as the value of the "fake" half-map. These are then useul for map masking, as the correlation between the "fake" half-maps and the original maps give nice masks according to Rangana
[in] | map | A Reference Pointer to the map which should be blurred/sharpened. |
[in] | blurredMap | A Reference Pointer to the variable which will store the modified map. |
[in] | xDimS | The number of indices along the x axis of the map. |
[in] | yDimS | The number of indices along the y axis of the map. |
[in] | zDimS | The number of indices along the z axis of the map. |
[in] | fakeMapKernel | The amount of neighbours in any direction whose average is to be used to get the current point. |
Definition at line 1680 of file ProSHADE_mapManip.cpp.
proshade_single ProSHADE_internal_mapManip::getIndicesFromAngstroms | ( | proshade_unsign | xDim, |
proshade_unsign | yDim, | ||
proshade_unsign | zDim, | ||
proshade_single | xAngs, | ||
proshade_single | yAngs, | ||
proshade_single | zAngs, | ||
proshade_single | dist | ||
) |
This function converts distance in Angstroms to distance in map indices.
This function finds out how many indices are required to cover a space of size "dist" in Angstroms. If you need this rounded in any way (ceil, floor, ...), just apply appropriate function to the output of this function.
[in] | xDim | The number of map indices along the x-axis. |
[in] | yDim | The number of map indices along the y-axis. |
[in] | zDim | The number of map indices along the z-axis. |
[in] | xAngs | The map size in Angstroms along the x-axis. |
[in] | yAngs | The map size in Angstroms along the y-axis. |
[in] | zAngs | The map size in Angstroms along the z-axis. |
[in] | dist | The distance in Angstroms to be converted |
Definition at line 1816 of file ProSHADE_mapManip.cpp.
void ProSHADE_internal_mapManip::getMaskFromBlurr | ( | proshade_double *& | blurMap, |
proshade_double *& | outMap, | ||
proshade_unsign | xDim, | ||
proshade_unsign | yDim, | ||
proshade_unsign | zDim, | ||
proshade_single | noIQRs | ||
) |
Function for computing mask from blurred map.
This function takes a blurred map and proceeds to compute the cut-off threshold for masking. It then applies this mask to the blurred map and applies the masking to the second input map. The assumption is that this second map is the map before blurring, as non-zero mask points will not be changed. Careful about this!!! It, however, does not output the mask.
[in] | blurMap | A Reference Pointer to the map which has been blurred for mask computation. |
[in] | outMap | A Reference Pointer to the map which will be masked - this map should be the map before blurring. |
[in] | xDim | The number of indices along the x axis of the map. |
[in] | yDim | The number of indices along the y axis of the map. |
[in] | zDim | The number of indices along the z axis of the map. |
[in] | noIQRs | The number of inter-quartile ranges from the median which should be used to compute the threshold for masking. |
Definition at line 1032 of file ProSHADE_mapManip.cpp.
void ProSHADE_internal_mapManip::getNonZeroBounds | ( | proshade_double * | map, |
proshade_signed | xDim, | ||
proshade_signed | yDim, | ||
proshade_signed | zDim, | ||
proshade_signed *& | ret | ||
) |
Function for finding the map boundaries enclosing positive only values.
This function takes a map and searches for minimum and maximum map indices in all three dimensions, which enclose all the non-zero map values.
[in] | map | A pointer to the map for which the bounds are to be found. |
[in] | xDim | The number of indices along the x axis of the map. |
[in] | yDim | The number of indices along the y axis of the map. |
[in] | zDim | The number of indices along the z axis of the map. |
[in] | ret | A proshade_unsign pointer to array of 6 to which the results will be saved (0 = minX; 1 = maxX; 2 = minY; 3 = maxY; 4 - minZ; 5 = maxZ). |
Definition at line 1082 of file ProSHADE_mapManip.cpp.
void ProSHADE_internal_mapManip::moveMapByFourier | ( | proshade_double *& | map, |
proshade_single | xMov, | ||
proshade_single | yMov, | ||
proshade_single | zMov, | ||
proshade_single | xAngs, | ||
proshade_single | yAngs, | ||
proshade_single | zAngs, | ||
proshade_signed | xDim, | ||
proshade_signed | yDim, | ||
proshade_signed | zDim | ||
) |
Function for moving map back to original PDB location by using Fourier transformation.
This function translates the map by changing the phase information of the map Fourier transform and then doing the inverse Fourier. This allows for movement smaller than one index distance, but should not be used over long distances (could move out of boundaries and meet pariodicity problem such as map from different cell now being moved into this cell).
[in] | map | A Reference Pointer to the map which should be shifted. |
[in] | xMov | The NEGATIVE value by how many angstroms should the x axis be moved. |
[in] | yMov | The NEGATIVE value by how many angstroms should the y axis be moved. |
[in] | zMov | The NEGATIVE value by how many angstroms should the z axis be moved. |
[in] | xAngs | How many angstroms are there along the x dimension. |
[in] | yAngs | How many angstroms are there along the y dimensionzAngs |
[in] | zAngs | How many angstroms are there along the z dimension. |
[in] | xDim | How many indices are there along the x dimension. |
[in] | yDim | How many indices are there along the y dimension. |
[in] | zDim | How many indices are there along the z dimension. |
Definition at line 811 of file ProSHADE_mapManip.cpp.
void ProSHADE_internal_mapManip::moveMapByIndices | ( | proshade_single * | xMov, |
proshade_single * | yMov, | ||
proshade_single * | zMov, | ||
proshade_single | xAngs, | ||
proshade_single | yAngs, | ||
proshade_single | zAngs, | ||
proshade_signed * | xFrom, | ||
proshade_signed * | xTo, | ||
proshade_signed * | yFrom, | ||
proshade_signed * | yTo, | ||
proshade_signed * | zFrom, | ||
proshade_signed * | zTo, | ||
proshade_signed * | xOrigin, | ||
proshade_signed * | yOrigin, | ||
proshade_signed * | zOrigin | ||
) |
Function for moving map back to original PDB location by changing the indices.
This function translates the map by changing the to and from index values so that the location of the map will be the same as the location of the original PDB file. This most likely cannot be done exactly as it allows only movement by increments of the sampling rate, but it is a decent start.
[in] | xMov | Pointer to the NEGATIVE value by how many angstroms should the x axis be moved. |
[in] | yMov | Pointer to the NEGATIVE value by how many angstroms should the y axis be moved. |
[in] | zMov | Pointer to the NEGATIVE value by how many angstroms should the z axis be moved. |
[in] | xAngs | How many angstroms are there along the x dimension. |
[in] | yAngs | How many angstroms are there along the y dimension. |
[in] | zAngs | How many angstroms are there along the z dimension. |
[in] | xFrom | The initial index of the x dimension of the map. |
[in] | xTo | The terminal index of the x dimension of the map. |
[in] | yFrom | The initial index of the y dimension of the map. |
[in] | yTo | The terminal index of the y dimension of the map. |
[in] | zFrom | The initial index of the z dimension of the map. |
[in] | zTo | The terminal index of the z dimension of the map. |
[in] | xOrigin | The first value of the x axis index. |
[in] | yOrigin | The first value of the y axis index. |
[in] | zOrigin | The first value of the z axis index. |
Definition at line 763 of file ProSHADE_mapManip.cpp.
void ProSHADE_internal_mapManip::movePDBForMapCalc | ( | gemmi::Structure * | pdbFile, |
proshade_single | xMov, | ||
proshade_single | yMov, | ||
proshade_single | zMov, | ||
bool | firstModel | ||
) |
Function for moving co-ordinate atoms to better suit theoretical map computation.
This function translates all atoms by a given x, y and z distances. This is required as theoretical map computation can only output map cells starting from 0, 0, 0 and therefore to avoid density being in corners for PDB atoms not located in posivite axes, the atoms need to be moved. This effect should be reversed later.
[in] | pdbFile | A pointer to the gemmi::Structure object as read in from the input file. |
[in] | xMov | How many angstroms should the atoms be moved along the x axis. |
[in] | yMov | How many angstroms should the atoms be moved along the y axis. |
[in] | zMov | How many angstroms should the atoms be moved along the z axis. |
[in] | firstModel | Should only the first, or all models be used? |
Definition at line 573 of file ProSHADE_mapManip.cpp.
proshade_signed ProSHADE_internal_mapManip::myRound | ( | proshade_double | x | ) |
Calls the appropriate version of round function depending on compiler version.
[in] | x | A decimal point number to be rounded. |
[out] | X | The rounded number. |
Definition at line 31 of file ProSHADE_mapManip.cpp.
proshade_signed ProSHADE_internal_mapManip::myRound | ( | proshade_single | x | ) |
Calls the appropriate version of round function depending on compiler version.
[in] | x | A decimal point number to be rounded. |
[out] | X | The rounded number. |
Definition at line 46 of file ProSHADE_mapManip.cpp.
void ProSHADE_internal_mapManip::releaseResolutionFourierMemory | ( | fftw_complex *& | origMap, |
fftw_complex *& | fCoeffs, | ||
fftw_complex *& | newFCoeffs, | ||
fftw_complex *& | newMap, | ||
fftw_plan & | planForwardFourier, | ||
fftw_plan & | planBackwardRescaledFourier | ||
) |
This function releases the memory required by the Fourier resampling.
This function simply deletes all the memory allocated by the allocateResolutionFourierMemory() function.
[in] | origMap | A Reference pointer to an array where the original map data were stored. |
[in] | fCoeffs | A Reference pointer to an array where the original map Fourier coefficients data were stored. |
[in] | newFCoeffs | A Reference pointer to an array where the re-sampled map Fourier coefficients data were stored. |
[in] | newMap | A Reference pointer to an array where the re-sampled map data were stored. |
[in] | planForwardFourier | FFTW_plan which computed the original map to Fourier coefficients transform. |
[in] | planBackwardRescaledFourier | FFTW_plan which computed the re-sampled Fourier coefficients to re-sampled map transform. |
Definition at line 1530 of file ProSHADE_mapManip.cpp.
void ProSHADE_internal_mapManip::removeMapPhase | ( | fftw_complex *& | mapCoeffs, |
proshade_unsign | xDim, | ||
proshade_unsign | yDim, | ||
proshade_unsign | zDim | ||
) |
This function removes the phase from reciprocal (frequency) map.
This function takes an already FFTW-ed map and its dimensions as the input and proceeds to remove the phase from the map. It writes over the map and does not release any memory - it is the role of the calling function to deal with both these features.
[in] | mapCoeffs | A Reference Pointer to the frequency map, from which phase is to be removed. |
[in] | xDim | The number of indices along the x-axis of the input map. |
[in] | yDim | The number of indices along the y-axis of the input map. |
[in] | zDim | The number of indices along the z-axis of the input map. |
Definition at line 1631 of file ProSHADE_mapManip.cpp.
void ProSHADE_internal_mapManip::removeWaters | ( | gemmi::Structure * | pdbFile, |
bool | firstModel | ||
) |
This function removed all waters from PDB input file.
This function does a quick read-through the PDB file and deletes all water molecules from the PDB file. This is important as comparing two similar structures, one with hundreds of waters (with high occupancy in some cases) and one without will lead to the structures being different in ProSHADE's eyes.
[in] | pdbFile | A pointer to gemmi::Structure object read in from the input file. |
[in] | firstModel | Should only the first, or all models be used? |
Definition at line 504 of file ProSHADE_mapManip.cpp.
void ProSHADE_internal_mapManip::reSampleMapToResolutionFourier | ( | proshade_double *& | map, |
proshade_single | resolution, | ||
proshade_unsign | xDimS, | ||
proshade_unsign | yDimS, | ||
proshade_unsign | zDimS, | ||
proshade_single | xAngs, | ||
proshade_single | yAngs, | ||
proshade_single | zAngs, | ||
proshade_single *& | corrs | ||
) |
This function re-samples a map to conform to given resolution using Fourier.
This function re-samples the internal map to a given resolutution by removing or zero-padding the Fourier (reciprocal space) coefficients and computing the inverse Fourier transform. This is the default option for map re-sampling, should it be required by the user.
[in] | map | A Reference Pointer to the map for which the bounds are to be found. |
[in] | resolution | The required resolution value. |
[in] | xDimS | The number of indices along the x axis of the map. |
[in] | yDimS | The number of indices along the y axis of the map. |
[in] | zDimS | The number of indices along the z axis of the map. |
[in] | xAngs | The size of the x dimension of the map in angstroms. |
[in] | yAngs | The size of the y dimension of the map in angstroms. |
[in] | zAngs | The size of the z dimension of the map in angstroms. |
[in] | corrs | Pointer reference to proshade_single array of 6 values with the following meaning: 0 = xAdd; 1 = yAdd; 2 = zAdd; 3 = newXAng; 4 = newYAng; 5 = newZAng |
Definition at line 1376 of file ProSHADE_mapManip.cpp.
void ProSHADE_internal_mapManip::reSampleMapToResolutionTrilinear | ( | proshade_double *& | map, |
proshade_single | resolution, | ||
proshade_unsign | xDimS, | ||
proshade_unsign | yDimS, | ||
proshade_unsign | zDimS, | ||
proshade_single | xAngs, | ||
proshade_single | yAngs, | ||
proshade_single | zAngs, | ||
proshade_single *& | corrs | ||
) |
This function re-samples a map to conform to given resolution using tri-linear interpolation.
This function takes a map and resolution value and it proceeds to create a new map, which has sampling resolution/2 and is large enough to contain the original map. It then proceeds to interpolate the new map values from the old map values, re-writing the old map once interpolation is done.
[in] | map | A Reference Pointer to the map for which the bounds are to be found. |
[in] | resolution | The required resolution value. |
[in] | xDimS | The number of indices along the x axis of the map. |
[in] | yDimS | The number of indices along the y axis of the map. |
[in] | zDimS | The number of indices along the z axis of the map. |
[in] | xAngs | The size of the x dimension of the map in angstroms. |
[in] | yAngs | The size of the y dimension of the map in angstroms. |
[in] | zAngs | The size of the z dimension of the map in angstroms. |
[in] | corrs | Pointer reference to proshade_single array of 6 values with the following meaning: 0 = xAdd; 1 = yAdd; 2 = zAdd; 3 = newXAng; 4 = newYAng; 5 = newZAng |
Definition at line 1175 of file ProSHADE_mapManip.cpp.
void ProSHADE_internal_mapManip::rotatePDBCoordinates | ( | gemmi::Structure * | pdbFile, |
proshade_double | euA, | ||
proshade_double | euB, | ||
proshade_double | euG, | ||
proshade_double | xCom, | ||
proshade_double | yCom, | ||
proshade_double | zCom, | ||
bool | firstModel | ||
) |
Function for rotating the PDB file co-ordinates by Euler angles.
This function takes the three Euler angles and a pointer to a gemmi::Structure and it then proceeds to compute the rotation matrix from the Euler angles. This matrix is then applied to the co-ordinates in the gemmi::Structure in a way so that the rotation is done over the centre of the co-ordinates, but the co-ordinate positions stay unchanged.
[in] | pdbFile | Pointer to a gemmi::Structure object which will have its co-ordinates rotated. |
[in] | euA | The Euler angle alpha by which the co-ordinates should be rotated. |
[in] | euB | The Euler angle beta by which the co-ordinates should be rotated. |
[in] | euG | The Euler angle gamma by which the co-ordinates should be rotated. |
[in] | xCom | The x-axis position around which the rotation should be applied. |
[in] | yCom | The y-axis position around which the rotation should be applied. |
[in] | zCom | The z-axis position around which the rotation should be applied. |
[in] | firstModel | Should only the first, or all models be used? |
Definition at line 296 of file ProSHADE_mapManip.cpp.
void ProSHADE_internal_mapManip::translatePDBCoordinates | ( | gemmi::Structure * | pdbFile, |
proshade_double | transX, | ||
proshade_double | transY, | ||
proshade_double | transZ, | ||
bool | firstModel | ||
) |
Function for translating the PDB file co-ordinates by given distances in Angstroms.
This function simply iterates through the given structure object and adds the required shift to all atoms in the first model along the three axes.
[in] | pdbFile | Pointer to a gemmi::Structure object whose co-ordinates are to be translated. |
[in] | transX | The |
[in] | transY | The |
[in] | transZ | The |
[in] | firstModel | Should only the first, or all models be used? |
Definition at line 381 of file ProSHADE_mapManip.cpp.