ProSHADE  0.7.6.6 (JUL 2022)
Protein Shape Detection
ProSHADE_mapManip.hpp
Go to the documentation of this file.
1 
22 //==================================================== ProSHADE
23 #include "ProSHADE_io.hpp"
24 
25 //==================================================== Overinclusion protection
26 #ifndef PROSHADE_MAPMANIP
27 #define PROSHADE_MAPMANIP
28 
29 //==================================================== ProSHADE_internal_mapManip namespace
38 {
39  proshade_signed myRound ( proshade_double x );
40  proshade_signed myRound ( proshade_single x );
41  void determinePDBRanges ( gemmi::Structure pdbFile, proshade_single* xFrom, proshade_single* xTo, proshade_single* yFrom,
42  proshade_single* yTo, proshade_single* zFrom, proshade_single* zTo, bool firstModel );
43  void findPDBCOMValues ( gemmi::Structure* pdbFile, proshade_double *xCom, proshade_double *yCom, proshade_double *zCom, bool firstModel );
44  void findMAPCOMValues ( proshade_double* map, proshade_double *xCom, proshade_double *yCom, proshade_double *zCom,
45  proshade_single xAngs, proshade_single yAngs, proshade_single zAngs, proshade_signed xFrom,
46  proshade_signed xTo, proshade_signed yFrom, proshade_signed yTo, proshade_signed zFrom,
47  proshade_signed zTo, bool removeNegDens );
48  void rotatePDBCoordinates ( gemmi::Structure *pdbFile, proshade_double euA, proshade_double euB, proshade_double euG, proshade_double xCom,
49  proshade_double yCom, proshade_double zCom, bool firstModel );
50  void translatePDBCoordinates ( gemmi::Structure *pdbFile, proshade_double transX, proshade_double transY, proshade_double transZ, bool firstModel );
51  void changePDBBFactors ( gemmi::Structure *pdbFile, proshade_double newBFactorValue, bool firstModel );
52  void removeWaters ( gemmi::Structure *pdbFile, bool firstModel );
53  void movePDBForMapCalc ( gemmi::Structure *pdbFile, proshade_single xMov, proshade_single yMov, proshade_single zMov, bool firstModel );
54  void generateMapFromPDB ( gemmi::Structure pdbFile, proshade_double*& map, proshade_single requestedResolution,
55  proshade_single xCell, proshade_single yCell, proshade_single zCell, proshade_signed* xTo,
56  proshade_signed* yTo, proshade_signed* zTo, bool forceP1, bool firstModel );
57  void moveMapByIndices ( proshade_single* xMov, proshade_single* yMov, proshade_single* zMov, proshade_single xAngs, proshade_single yAngs,
58  proshade_single zAngs, proshade_signed* xFrom, proshade_signed* xTo, proshade_signed* yFrom, proshade_signed* yTo,
59  proshade_signed* zFrom, proshade_signed* zTo, proshade_signed* xOrigin, proshade_signed* yOrigin,
60  proshade_signed* zOrigin );
61  void moveMapByFourier ( proshade_double*& map, proshade_single xMov, proshade_single yMov, proshade_single zMov, proshade_single xAngs,
62  proshade_single yAngs, proshade_single zAngs, proshade_signed xDim, proshade_signed yDim, proshade_signed zDim );
63  void moveMapByFourierInReci ( proshade_complex*& coeffs, proshade_double*& weights, proshade_single xMov, proshade_single yMov, proshade_single zMov, proshade_single xAngs,
64  proshade_single yAngs, proshade_single zAngs, proshade_signed xDim, proshade_signed yDim, proshade_signed zDim );
65  void blurSharpenMap ( proshade_double*& map, proshade_double*& maskedMap, proshade_unsign xDimS, proshade_unsign yDimS,
66  proshade_unsign zDimS, proshade_single xAngs, proshade_single yAngs, proshade_single zAngs, proshade_single blurringFactor );
67  void getMaskFromBlurr ( proshade_double*& blurMap, proshade_double*& outMap, proshade_unsign xDimS, proshade_unsign yDimS,
68  proshade_unsign zDimS, proshade_single noIQRs );
69  void getNonZeroBounds ( proshade_double* map, proshade_signed xDim, proshade_signed yDim, proshade_signed zDim, proshade_signed*& ret );
70  void addExtraBoundSpace ( proshade_unsign xDim, proshade_unsign yDim, proshade_unsign zDim, proshade_single xAngs, proshade_single yAngs,
71  proshade_single zAngs, proshade_signed*& bounds, proshade_single extraSpace );
72  void reSampleMapToResolutionTrilinear ( proshade_double*& map, proshade_single resolution, proshade_unsign xDimS, proshade_unsign yDimS,
73  proshade_unsign zDimS, proshade_single xAngs, proshade_single yAngs, proshade_single zAngs,
74  proshade_single*& corrs );
75  void reSampleMapToResolutionFourier ( proshade_double*& map, proshade_single resolution, proshade_unsign xDimS, proshade_unsign yDimS,
76  proshade_unsign zDimS, proshade_single xAngs, proshade_single yAngs, proshade_single zAngs,
77  proshade_single*& corrs );
78  void allocateResolutionFourierMemory ( fftw_complex*& origMap, fftw_complex*& fCoeffs, fftw_complex*& newFCoeffs, fftw_complex*& newMap, fftw_plan& planForwardFourier,
79  fftw_plan& planBackwardRescaledFourier, proshade_unsign xDimOld, proshade_unsign yDimOld, proshade_unsign zDimOld,
80  proshade_unsign xDimNew, proshade_unsign yDimNew, proshade_unsign zDimNew );
81  void releaseResolutionFourierMemory ( fftw_complex*& origMap, fftw_complex*& fCoeffs, fftw_complex*& newFCoeffs, fftw_complex*& newMap, fftw_plan& planForwardFourier,
82  fftw_plan& planBackwardRescaledFourier );
83  void changeFourierOrder ( fftw_complex*& fCoeffs, proshade_signed xDim, proshade_signed yDim, proshade_signed zDim, bool negativeFirst );
84  void removeMapPhase ( fftw_complex*& mapCoeffs, proshade_unsign xDim, proshade_unsign yDim, proshade_unsign zDim );
85  void getFakeHalfMap ( proshade_double*& map, proshade_double*& fakeHalfMap, proshade_unsign xDimS, proshade_unsign yDimS,
86  proshade_unsign zDimS, proshade_signed fakeMapKernel );
87  void getCorrelationMapMask ( proshade_double*& map, proshade_double*& fakeHalfMap, proshade_double*& correlationMask,
88  proshade_unsign xDimS, proshade_unsign yDimS,
89  proshade_unsign zDimS, proshade_signed corrMaskKernel );
90  proshade_single getIndicesFromAngstroms ( proshade_unsign xDim, proshade_unsign yDim, proshade_unsign zDim,
91  proshade_single xAngs, proshade_single yAngs, proshade_single zAngs,
92  proshade_single dist );
93  void connectMaskBlobs ( proshade_double*& mask, proshade_signed xDim, proshade_signed yDim, proshade_signed zDim,
94  proshade_single xAngs, proshade_single yAngs,
95  proshade_single zAngs, proshade_single maskThres );
96  void beautifyBoundaries ( proshade_signed*& bounds, proshade_unsign xDim, proshade_unsign yDim, proshade_unsign zDim, proshade_signed boundsDiffThres );
97  proshade_signed betterClosePrimeFactors ( proshade_signed fromRange, proshade_signed toRange );
98  void distributeSpaceToBoundaries ( proshade_signed& minBound, proshade_signed& maxBound, proshade_signed oldBoundRange, proshade_signed newBoundRange );
99  void copyMapByBounds ( proshade_signed xFrom, proshade_signed xTo, proshade_signed yFrom, proshade_signed yTo, proshade_signed zFrom,
100  proshade_signed zTo, proshade_signed origXFrom, proshade_signed origYFrom, proshade_signed origZFrom,
101  proshade_unsign yDimIndices, proshade_unsign zDimIndices, proshade_unsign origXDimIndices,
102  proshade_unsign origYDimIndices, proshade_unsign origZDimIndices, proshade_double*& newMap,
103  proshade_double* origMap );
104 }
105 
106 #endif
ProSHADE_internal_mapManip::addExtraBoundSpace
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.
Definition: ProSHADE_mapManip.cpp:1193
ProSHADE_internal_mapManip::determinePDBRanges
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.
Definition: ProSHADE_mapManip.cpp:72
ProSHADE_internal_mapManip::changePDBBFactors
void changePDBBFactors(gemmi::Structure *pdbFile, proshade_double newBFactorValue, bool firstModel)
Function for changing the PDB B-factor values to a specific single value.
Definition: ProSHADE_mapManip.cpp:458
ProSHADE_internal_mapManip::distributeSpaceToBoundaries
void distributeSpaceToBoundaries(proshade_signed &minBound, proshade_signed &maxBound, proshade_signed oldBoundRange, proshade_signed newBoundRange)
Function for adding space to boundaries within the map confines.
Definition: ProSHADE_mapManip.cpp:2092
ProSHADE_internal_mapManip::getMaskFromBlurr
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.
Definition: ProSHADE_mapManip.cpp:1086
ProSHADE_internal_mapManip::blurSharpenMap
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.
Definition: ProSHADE_mapManip.cpp:986
ProSHADE_internal_mapManip::moveMapByFourier
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.
Definition: ProSHADE_mapManip.cpp:825
ProSHADE_internal_mapManip::movePDBForMapCalc
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.
Definition: ProSHADE_mapManip.cpp:587
ProSHADE_internal_mapManip::reSampleMapToResolutionTrilinear
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.
Definition: ProSHADE_mapManip.cpp:1229
ProSHADE_io.hpp
This header file declares all the functions required for low level file format access.
ProSHADE_internal_mapManip::getCorrelationMapMask
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.
Definition: ProSHADE_mapManip.cpp:1805
ProSHADE_internal_mapManip::changeFourierOrder
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...
Definition: ProSHADE_mapManip.cpp:1618
ProSHADE_internal_mapManip::releaseResolutionFourierMemory
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.
Definition: ProSHADE_mapManip.cpp:1589
ProSHADE_internal_mapManip::removeWaters
void removeWaters(gemmi::Structure *pdbFile, bool firstModel)
This function removed all waters from PDB input file.
Definition: ProSHADE_mapManip.cpp:518
ProSHADE_internal_mapManip::getIndicesFromAngstroms
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.
Definition: ProSHADE_mapManip.cpp:1876
ProSHADE_internal_mapManip::getNonZeroBounds
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.
Definition: ProSHADE_mapManip.cpp:1136
ProSHADE_internal_mapManip::findMAPCOMValues
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, bool removeNegDens)
This function finds the Centre of Mass for a map.
Definition: ProSHADE_mapManip.cpp:243
ProSHADE_internal_mapManip
This namespace contains the internal functions for manipulating maps already present in the internal ...
Definition: ProSHADE_io.cpp:27
ProSHADE_internal_mapManip::removeMapPhase
void removeMapPhase(fftw_complex *&mapCoeffs, proshade_unsign xDim, proshade_unsign yDim, proshade_unsign zDim)
This function removes the phase from reciprocal (frequency) map.
Definition: ProSHADE_mapManip.cpp:1691
ProSHADE_internal_mapManip::allocateResolutionFourierMemory
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.
Definition: ProSHADE_mapManip.cpp:1555
ProSHADE_internal_mapManip::betterClosePrimeFactors
proshade_signed betterClosePrimeFactors(proshade_signed fromRange, proshade_signed toRange)
Function for finding close numbers with better prime factors.
Definition: ProSHADE_mapManip.cpp:2046
ProSHADE_internal_mapManip::myRound
proshade_signed myRound(proshade_double x)
Calls the appropriate version of round function depending on compiler version.
Definition: ProSHADE_mapManip.cpp:31
ProSHADE_internal_mapManip::findPDBCOMValues
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.
Definition: ProSHADE_mapManip.cpp:157
ProSHADE_internal_mapManip::beautifyBoundaries
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.
Definition: ProSHADE_mapManip.cpp:1976
ProSHADE_internal_mapManip::rotatePDBCoordinates
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.
Definition: ProSHADE_mapManip.cpp:310
ProSHADE_internal_mapManip::generateMapFromPDB
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.
Definition: ProSHADE_mapManip.cpp:658
ProSHADE_internal_mapManip::translatePDBCoordinates
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.
Definition: ProSHADE_mapManip.cpp:395
ProSHADE_internal_mapManip::getFakeHalfMap
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.
Definition: ProSHADE_mapManip.cpp:1740
ProSHADE_internal_mapManip::reSampleMapToResolutionFourier
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.
Definition: ProSHADE_mapManip.cpp:1439
ProSHADE_internal_mapManip::copyMapByBounds
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.
Definition: ProSHADE_mapManip.cpp:2142
ProSHADE_internal_mapManip::connectMaskBlobs
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.
Definition: ProSHADE_mapManip.cpp:1899
ProSHADE_internal_mapManip::moveMapByFourierInReci
void moveMapByFourierInReci(proshade_complex *&coeffs, proshade_double *&weights, 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 using Fourier coefficients in the reciprocal space only.
Definition: ProSHADE_mapManip.cpp:910
ProSHADE_internal_mapManip::moveMapByIndices
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.
Definition: ProSHADE_mapManip.cpp:777