ProSHADE  0.6.6 (DEC 2018)
Protein Shape Descriptors and Symmetry Detection
ProSHADE_internal.h
Go to the documentation of this file.
1 
19 //============================================ ProSHADE
20 #include "ProSHADE_version.h"
21 #include "ProSHADE_rvapi.h"
22 
23 //============================================ Include once
24 #ifndef __PROSHADE_INTERNAL_LIBRARY__
25 #define __PROSHADE_INTERNAL_LIBRARY__
26 
33 namespace ProSHADE_internal
34 {
35  //======================================== Forward declaration
36  unsigned int checkFileType ( std::string fileName );
37 
38  //======================================== Forward declaration
39  class ProSHADE_comparePairwise;
40 
49  {
50  //==================================== Friends
51  friend class ProSHADE_comparePairwise;
52 
53  private:
54  //==================================== Settings regarding resolutions
55  double mapResolution;
56  std::vector<unsigned int> bandwidth;
57  std::vector<unsigned int> glIntegOrder;
58  std::vector<unsigned int> theta;
59  std::vector<unsigned int> phi;
60 
61  //==================================== Settings regarding B factors
62  double bFactorValue;
63  double bFactorChange;
64 
65  //==================================== Setting regarding maps and removing noise
66  double noIQRsFromMap;
67 
68  //==================================== Settings regarding concentric shells
69  double shellSpacing;
70  unsigned int manualShells;
71 
72  //==================================== Settings regarding map with phases
73  bool useCOM;
74  bool firstLineCOM;
75 
76  //==================================== Settings regarding space around the structure in lattice
77  double extraSpace;
78 
79  //==================================== Settings regarding weighting the distances
80  double alpha;
81  double mPower;
82 
83  //==================================== Settings regarding bands to ignore
84  std::vector<int> ignoreLs;
85 
86  //==================================== Settings regarding the structures to use
87  std::vector <std::string> structFiles;
88 
89  //==================================== Variables holding the results
90  std::vector< std::array<double,8> > rfPeaks;
91  std::vector< std::vector<std::array<double,5> > > tetrSymm;
92  double tetrSymmPeakAvg;
93  std::vector< std::array<double,5> > tetrElems;
94  std::vector< std::vector<std::array<double,5> > > octaSymm;
95  double octaSymmPeakAvg;
96  std::vector< std::array<double,5> > octaElems;
97  std::vector< std::vector<std::array<double,5> > > icosSymm;
98  double icosSymmPeakAvg;
99  std::vector< std::array<double,5> > icosElems;
100 
101  //==================================== Settings regarding peak finding
102  double peakHeightNoIQRs;
103  double peakDistanceForReal;
104  int peakSurroundingPoints;
105 
106  //==================================== Settings regarding axis and angle error tolerance
107  double aaErrorTolerance;
108  double symGapTolerance;
109 
110  //==================================== Are full results needed?
111  bool printFull;
112 
113  public:
114  //==================================== Public variables holding the results
115  std::vector< std::array<double,5> > cnSymm;
116  std::vector< std::array<double,5> > cnSymmClear;
117  std::vector< std::vector<std::array<double,6> > > dnSymm;
118  std::vector< std::vector<std::array<double,6> > > dnSymmClear;
119  std::vector< std::array<double,5> > tetrAxes;
120  std::vector< std::array<double,5> > octaAxes;
121  std::vector< std::array<double,5> > icosAxes;
123 
124  //==================================== Public functions
126  ProSHADE_symmetry ( );
127  void printResultsClear ( int verbose );
129  void printResultsRequest ( std::string symmetryType,
130  unsigned int symmetryFold,
131  int verbose );
132  void printResultsRequestHTML ( std::string symmetryType,
133  unsigned int symmetryFold,
134  ProSHADE::ProSHADE_settings* settings );
135  std::vector< std::array<double,8> > getRotFnPeaks ( );
136  std::vector< std::array<double,5> > getCSymmetries ( );
137  std::vector< std::vector<std::array<double,6> > > getDSymmetries ( );
138  std::vector< std::array<double,5> > generateTetrAxes ( ProSHADE_comparePairwise* cmpObj,
139  ProSHADE::ProSHADE_settings* settings,
140  std::vector< std::array<double,5> > tetrSymm,
141  std::vector< std::array<double,5> > allCs,
142  double axisErrorTolerance = 0.1,
143  int verbose = 0 );
144  std::vector< std::array<double,5> > generateTetrElements ( std::vector< std::array<double,5> > symmAxes,
145  ProSHADE::ProSHADE_settings* settings,
146  int verbose = 0 );
147  std::vector< std::array<double,5> > generateOctaAxes ( ProSHADE_comparePairwise* cmpObj,
148  ProSHADE::ProSHADE_settings* settings,
149  std::vector< std::array<double,5> > octaSymm,
150  std::vector< std::array<double,5> > allCs,
151  double axisErrorTolerance = 0.1,
152  int verbose = 0 );
153  std::vector< std::array<double,5> > generateOctaElements ( std::vector< std::array<double,5> > symmAxes,
154  ProSHADE::ProSHADE_settings* settings,
155  int verbose = 0 );
156  std::vector< std::array<double,5> > generateIcosAxes ( ProSHADE_comparePairwise* cmpObj,
157  ProSHADE::ProSHADE_settings* settings,
158  std::vector< std::array<double,5> > icosSymm,
159  std::vector< std::array<double,5> > allCs,
160  double axisErrorTolerance = 0.1,
161  int verbose = 0 );
162  std::vector< std::array<double,5> > generateIcosElements ( std::vector< std::array<double,5> > symmAxes,
163  ProSHADE::ProSHADE_settings* settings,
164  int verbose = 0 );
165  };
166 
175  {
176  //==================================== Friends
177  friend class ProSHADE_comparePairwise;
178  friend class ProSHADE_compareOneAgainstAll;
179  friend class ProSHADE_symmetry;
180  friend class ProSHADE_distances;
181  friend class ProSHADE_mapFeatures;
182 
183  private:
184  //==================================== Variables regarding the reading in of the data
185  std::string _inputFileName;
186  bool _fromPDB;
187 
188  //==================================== Variables regarding the shells
189  double _shellSpacing;
190  double _maxExtraCellularSpace;
191  int _exCeSpDiffX;
192  int _exCeSpDiffY;
193  int _exCeSpDiffZ;
194  double _xRange;
195  double _yRange;
196  double _zRange;
197  double _xSamplingRate;
198  double _ySamplingRate;
199  double _zSamplingRate;
200  float _xFrom;
201  float _yFrom;
202  float _zFrom;
203  float _xTo;
204  float _yTo;
205  float _zTo;
206  std::vector<double> _shellPlacement;
207 
208  //==================================== Variables regarding the density map
209  double _mapResolution;
210  unsigned int _maxMapU;
211  unsigned int _maxMapV;
212  unsigned int _maxMapW;
213  float *_densityMapMap;
214  bool _densityMapComputed;
215  double _mapMean;
216  double _mapSdev;
217  int _preCBSU;
218  int _preCBSV;
219  int _preCBSW;
220 
221  //==================================== Variables regarding the phase removal from map
222  double _fourierCoeffPower;
223  double _bFactorChange;
224  double *_densityMapCor;
225  double _maxMapRange;
226  std::array<double,3> *_densityMapCorCoords;
227  bool _phaseRemoved;
228  bool _usePhase;
229  bool _keepOrRemove;
230 
231  //==================================== Variables regarding the sphere mapping of phaseless data
232  double _thetaAngle;
233  double _phiAngle;
234  unsigned int _noShellsWithData;
235  double **_shellMappedData;
236  int _xCorrection;
237  int _yCorrection;
238  int _zCorrection;
239  double _xCorrErr;
240  double _yCorrErr;
241  double _zCorrErr;
242  bool _sphereMapped;
243  bool _firstLineCOM;
244 
245  //==================================== Variables regarding the spherical harmonics decomposition
246  unsigned int _bandwidthLimit;
247  unsigned int _oneDimmension;
248  double **_realSHCoeffs;
249  double **_imagSHCoeffs;
250  double *_sphericalHarmonicsWeights;
251  double **_semiNaiveTable;
252  double *_semiNaiveTableSpace;
253  fftw_complex *_shWorkspace;
254  bool _sphericalCoefficientsComputed;
255  bool _wasBandwithGiven;
256  bool _wasThetaGiven;
257  bool _wasPhiGiven;
258  bool _wasGlInterGiven;
259 
260  //==================================== Variables regarding spherical harmonics inverse
261  double **_invRealData;
262  double **_invImagData;
263 
264  //==================================== Variables regarding the energy levels descriptor pre-calculation
265  double ***_rrpMatrices;
266  bool _rrpMatricesPrecomputed;
267 
268  public:
269  //==================================== Public functions
270  ProSHADE_data ( );
271  ProSHADE_data ( ProSHADE_data *copyFrom );
272  ProSHADE_data ( std::ifstream* dbFile, std::string fName, double volThreMin, double volThreMax, int verbose, ProSHADE::ProSHADE_settings* settings, bool saveWithAndWithout = false );
273  ~ProSHADE_data ( );
274 
275  void getDensityMapFromPDB ( std::string fileName, double *shellDistance, double resolution,
276  unsigned int *bandwidth, unsigned int *theta,
277  unsigned int *phi, unsigned int *glIntegOrder,
278  double *extraSpace, bool mapResDefault, ProSHADE::ProSHADE_settings* settings,
279  double Bfactor = 80.0, bool hpFirstLineCom = false, bool overlayDefaults = false );
280  void getDensityMapFromMAP ( std::string fileName, double *shellDistance, double resolution,
281  unsigned int *bandwidth, unsigned int *theta,
282  unsigned int *phi, unsigned int *glIntegOrder,
283  double *extraSpace, bool mapResDefault, bool rotDefaults, ProSHADE::ProSHADE_settings* settings,
284  bool overlayDefaults = false );
285  void maskMap ( int hlpU, int hlpV, int hlpW, double blurBy = 250.0, double maxMapIQR = 3.0, double extraMaskSpace = 3.0 );
286  void normaliseMap ( ProSHADE::ProSHADE_settings* settings );
287  std::array<double,4> getCOMandDist ( ProSHADE::ProSHADE_settings* settings );
288  bool removeIslands ( int hlpU, int hlpV, int hlpW, int verbose = 0, bool runAll = false );
289  void translateMap ( double xShift, double yShift, double zShift );
290  void shiftMap ( double xShift, double yShift, double zShift );
291  std::array<double,6> getDensityMapFromMAPRebox ( std::string fileName, double maxMapIQR, double extraCS, int verbose,
292  bool useCubicMaps, double maskBlurFactor, bool maskBlurFactorGiven,
293  ProSHADE::ProSHADE_settings* settings );
294  void getDensityMapFromMAPFeatures ( std::string fileName, double *minDensPreNorm, double *maxDensPreNorm,
295  double *minDensPostNorm, double *maxDensPostNorm, double *postNormMean,
296  double *postNormSdev, double *maskVolume, double* totalVolume, double *maskMean,
297  double *maskSdev, double *maskMin, double *maskMax, double *maskDensityRMS,
298  double *allDensityRMS, std::array<double,3> *origRanges, std::array<double,3> *origDims,
299  double maxMapIQR, double extraCS, int verbose, bool useCubicMaps, double maskBlurFactor, bool maskBlurFactorGiven,
300  bool reboxAtAll, ProSHADE::ProSHADE_settings* settings );
301  void removePhaseFromMap ( double alpha, double bFac, ProSHADE::ProSHADE_settings* settings );
302  void removePhaseFromMapOverlay ( double alpha, double bFac, ProSHADE::ProSHADE_settings* settings );
303  void keepPhaseInMap ( double alpha, double bFac,
304  unsigned int *bandwidth, unsigned int *theta,
305  unsigned int *phi, unsigned int *glIntegOrder, ProSHADE::ProSHADE_settings* settings,
306  bool useCom = true,
307  double maxMapIQR = 10.0, int verbose = 0,
308  bool clearMapData = true, bool rotDefaults = false, bool overlapDefaults = false,
309  double blurFactor = 500.0, bool maskBlurFactorGiven = false );
310  void mapPhaselessToSphere ( ProSHADE::ProSHADE_settings* settings, double theta, double phi, double shellSz, unsigned int manualShells = 0, bool keepInMemory = false, bool rotDefaults = false );
311  void getSphericalHarmonicsCoeffs ( unsigned int bandwidth, ProSHADE::ProSHADE_settings* settings );
312  void precomputeRotInvDescriptor ( ProSHADE::ProSHADE_settings* settings );
313  std::vector<ProSHADE_data*> fragmentMap ( ProSHADE::ProSHADE_settings* settings, bool userCOM );
314  std::vector< std::vector<int> > findMapIslands ( int hlpU, int hlpV, int hlpW, double *map, double threshold = 0.0 );
315  void writeMap ( std::string fileName,
316  double *map,
317  std::string axOrder = "xyz",
318  float xFrom = std::numeric_limits<float>::infinity(),
319  float yFrom = std::numeric_limits<float>::infinity(),
320  float zFrom = std::numeric_limits<float>::infinity(),
321  float xTo = std::numeric_limits<float>::infinity(),
322  float yTo = std::numeric_limits<float>::infinity(),
323  float zTo = std::numeric_limits<float>::infinity(),
324  float xRange = std::numeric_limits<float>::infinity(),
325  float yRange = std::numeric_limits<float>::infinity(),
326  float zRange = std::numeric_limits<float>::infinity() );
327  void writeMap ( std::string fileName,
328  float *map,
329  std::string axOrder = "xyz",
330  float xFrom = std::numeric_limits<float>::infinity(),
331  float yFrom = std::numeric_limits<float>::infinity(),
332  float zFrom = std::numeric_limits<float>::infinity(),
333  float xTo = std::numeric_limits<float>::infinity(),
334  float yTo = std::numeric_limits<float>::infinity(),
335  float zTo = std::numeric_limits<float>::infinity(),
336  float xRange = std::numeric_limits<float>::infinity(),
337  float yRange = std::numeric_limits<float>::infinity(),
338  float zRange = std::numeric_limits<float>::infinity() );
339  void writePDB ( std::string templateName, std::string outputName,
340  double rotEulA = 0.0,
341  double rotEulB = 0.0,
342  double rotEulG = 0.0,
343  double trsX = 0.0,
344  double trsY = 0.0,
345  double trsZ = 0.0 );
346  void deleteModel ( std::string modelPath );
347 
348  //==================================== Accessor functions
349  inline double getMapXRange ( void ) { return ( this->_xRange ); };
350  inline double getMapYRange ( void ) { return ( this->_yRange ); };
351  inline double getMapZRange ( void ) { return ( this->_zRange ); };
352  inline double* getMap ( void ) { return ( this->_densityMapCor ); };
353 
354  };
355 
364  {
365  private:
366  //==================================== Variables regarding sanity check
367  unsigned int _bandwidthLimit;
368  double _shellSpacing;
369  double _thetaAngle;
370  double _phiAngle;
371  double _matrixPowerWeight;
372  unsigned int _minShellsToUse;
373  std::vector<int> _lsToIgnore;
374  bool _bothRRPsPreComputed;
375  bool _keepOrRemove;
376 
377  //==================================== Variables regarding RotInv distance
378  double ***_obj1RRPs;
379  double ***_obj2RRPs;
380  double _distanceRotInv;
381  bool _rotInvComputed;
382 
383  //==================================== Variables regarding TrSigma distance pre-calculation
384  unsigned int _noShellsObj1;
385  unsigned int _noShellsObj2;
386  unsigned int _maxShellsToUse;
387  double **_obj1RealCoeffs;
388  double **_obj1ImagCoeffs;
389  double **_obj2RealCoeffs;
390  double **_obj2ImagCoeffs;
391  std::vector<double> _trSigmaWeights;
392  std::array<double,2> ***_trSigmaEMatrix;
393  bool _trSigmaPreComputed;
394 
395  //==================================== Variables regarding TrSigma distance computation
396  double _distanceTrSigma;
397  bool _trSigmaComputed;
398 
399  //==================================== Variables regarding the SO3 Transform
400  fftw_complex *_so3Coeffs;
401  fftw_complex *_so3InvCoeffs;
402  fftw_complex *_so3Workspace1;
403  fftw_complex *_so3Workspace2;
404  double *_so3Workspace3;
405  bool _so3InvMapComputed;
406 
407  //==================================== Variables regarding the Euler angles
408  std::array<double,3> _eulerAngles;
409  bool _eulerAnglesFound;
410 
411  //==================================== Variables regarding symmetry detection
412  double _peakHeightThr;
413  bool _CSymmsFound;
414  bool _DSymmsFound;
415 
416  //==================================== Variables regarding the Wigner d matrices computation
417  std::vector< std::vector< std::vector< std::array<double,2> > > > _wignerMatrices;
418  bool _wignerMatricesComputed;
419 
420  //==================================== Variables regarding Gauss-Legendre integration
421  unsigned int _glIntegrationOrder;
422  std::vector<double> _glAbscissas;
423  std::vector<double> _glWeights;
424 
425  //==================================== Variables regarding the coefficient rotation
426  bool _structureCoeffsRotated;
427 
428  //==================================== Internal functions
429  double gl20IntRR ( std::vector<double>* vals );
430  std::array<double,2> gl20IntCR ( std::vector< std::array<double,2> >* vals );
431  std::vector<double> getSingularValues ( unsigned int band, unsigned int dim, ProSHADE::ProSHADE_settings* settings );
432  std::vector< std::vector< std::vector<double> > > getSingularValuesUandVMatrix ( std::vector< std::vector<double> > mat, unsigned int dim, ProSHADE::ProSHADE_settings* settings );
433  double maxPeakHeightFromAngleAxis ( double X, double Y, double Z, double Angle, ProSHADE::ProSHADE_settings* settings );
434 
435  public:
436  //==================================== Public functions
437  ProSHADE_comparePairwise ( ProSHADE_data *cmpObj1, ProSHADE_data *cmpObj2, double mPower, std::vector<int> ignoreL, unsigned int order, ProSHADE::ProSHADE_settings* settings );
439 
440  void precomputeTrSigmaDescriptor ( );
441  void getSO3InverseMap ( ProSHADE::ProSHADE_settings* settings );
442  std::array<double,3> getEulerAngles ( ProSHADE::ProSHADE_settings* settings, double* correlation = nullptr );
443  void setEulerAngles ( double alpha, double beta, double gamma );
444  std::vector< std::array<double,8> > getSO3Peaks ( ProSHADE::ProSHADE_settings* settings, double noIQRs = 1.5, bool freeMem = true, int peakSize = 1, double realDist = 0.4, int verbose = 1 );
445  bool checkPeakExistence ( double alpha, double beta, double gamma, ProSHADE::ProSHADE_settings* settings, double passThreshold = 0.3 );
446  std::vector< std::array<double,5> > findCnSymmetry ( std::vector< std::array<double,8> > peaks,
447  ProSHADE::ProSHADE_settings* settings,
448  double axisErrorTolerance = 0.0,
449  bool freeMem = true,
450  double percAllowedToMiss = 0.33,
451  int verbose = 1 );
452  std::vector< std::array<double,5> > findCnSymmetryClear ( std::vector< std::array<double,5> > CnSymm, ProSHADE::ProSHADE_settings* settings, double maxGap = 0.2, bool *pf = nullptr );
453  std::vector< std::vector< std::array<double,6> > > findDnSymmetry ( std::vector< std::array<double,5> > CnSymm, double axisErrorTolerance = 0.1 );
454  std::vector< std::vector< std::array<double,6> > > findDnSymmetryClear ( std::vector< std::vector< std::array<double,6> > > DnSymm, ProSHADE::ProSHADE_settings* settings, double maxGap = 0.2, bool *pf = nullptr );
455  std::vector< std::vector< std::array<double,5> > > findTetrSymmetry ( std::vector< std::array<double,5> > CnSymm, double *tetrSymmPeakAvg, double axisErrorTolerance = 0.1 );
456  std::vector< std::vector< std::array<double,5> > > findOctaSymmetry ( std::vector< std::array<double,5> > CnSymm, double *octaSymmPeakAvg, double axisErrorTolerance = 0.1 );
457  std::vector< std::vector< std::array<double,5> > > findIcosSymmetry ( std::vector< std::array<double,5> > CnSymm, double *icosSymmPeakAvg, double axisErrorTolerance = 0.1 );
458  void generateWignerMatrices ( ProSHADE::ProSHADE_settings* settings );
459  double getRotCoeffDistance ( ProSHADE::ProSHADE_settings* settings );
460  double maxAvgPeakForSymmetry ( double X, double Y, double Z, double Angle, ProSHADE::ProSHADE_settings* settings );
461  void freeInvMap ( );
462  double getPeakThreshold ( );
463  void rotateStructure ( ProSHADE_data *cmpObj1, ProSHADE::ProSHADE_settings* settings,
464  std::string saveName, int verbose = 0, std::string axOrd = "xyz",
465  bool internalUse = false );
466  std::array<double,4> getTranslationFunctionMap ( ProSHADE_data *obj1, ProSHADE_data *obj2, double *ob2XMov = nullptr,
467  double *ob2YMov = nullptr, double *ob2ZMov = nullptr );
468  };
469 
478  {
479  //==================================== Friends
480  friend class ProSHADE_distances;
481 
482  private:
483  //==================================== Just save the input data and make the internal functions do the checks themselves
484  ProSHADE_data *one;
485  ProSHADE_data *two;
486  std::vector<ProSHADE_data*> *all;
487 
488  //==================================== These needs to be the same
489  std::vector<int> _lsToIgnore;
490  double _matrixPowerWeight;
491 
492  //==================================== Simple checks
493  bool _energyLevelsComputed;
494  bool _trSigmaPreComputed;
495  bool _trSigmaComputed;
496  bool _keepOrRemove;
497  bool _so3InvMapComputed;
498  bool _eulerAnglesFound;
499  bool _wignerMatricesComputed;
500  bool _fullDistComputed;
501 
502  //==================================== Result holders
503  std::vector<double> _distancesEnergyLevels;
504  std::vector<double> _distancesTraceSigma;
505  std::vector<double> _distancesFullRotation;
506 
507  //==================================== Internal value holders
508  std::vector< std::vector< std::vector< std::vector< std::array<double,2> > > > > _EMatrices;
509  std::vector<fftw_complex*> _so3InverseCoeffs;
510  std::vector< std::array<double,3> > _eulerAngles;
511  std::vector< std::vector< std::vector< std::vector< std::array<double,2> > > > > _wignerMatrices;
512  std::vector<unsigned int> _enLevelsDoNotFollow;
513  std::vector<unsigned int> _trSigmaDoNotFollow;
514 
515  //==================================== Internal functions
516  double glIntRR ( std::vector<double>* vals, double shellSep, std::vector<double>* glAbscissas, std::vector<double>* glWeights );
517  std::array<double,2> glIntCR ( std::vector< std::array<double,2> >* vals, double shellSep, std::vector<double> *glAbscissas, std::vector<double> *glWeights );
518  std::vector<double> getSingularValues ( unsigned int strNo, unsigned int band, unsigned int dim, ProSHADE::ProSHADE_settings* settings );
519 
520  public:
521  //==================================== Public functions
522  ProSHADE_compareOneAgainstAll ( ProSHADE_data *oneStr, std::vector<ProSHADE_data*> *allStrs,
523  std::vector<int> ignoreL, double matrixPowerWeight, int verbose );
524  ProSHADE_compareOneAgainstAll ( ProSHADE_data *oneStr, ProSHADE_data *twoStr, std::vector<ProSHADE_data*> *allStrs,
525  std::vector<int> ignoreL, double matrixPowerWeight, int verbose );
527 
528  std::vector<double> getEnergyLevelsDistance ( int verbose, ProSHADE::ProSHADE_settings* settings );
529  void precomputeTrSigmaDescriptor ( double shellSpacing, std::vector<unsigned int>* glIntegOrderVec, ProSHADE::ProSHADE_settings* settings );
530  std::vector<double> getTrSigmaDistance ( int verbose, ProSHADE::ProSHADE_settings* settings );
531  void getSO3InverseMap ( ProSHADE::ProSHADE_settings* settings );
532  std::vector< std::array<double,3> > getEulerAngles ( ProSHADE::ProSHADE_settings* settings );
533  void generateWignerMatrices ( ProSHADE::ProSHADE_settings* settings );
534  std::vector<double> getRotCoeffDistance ( int verbose, ProSHADE::ProSHADE_settings* settings );
535  void alignDensities ( ProSHADE::ProSHADE_settings* settings );
536  };
537 
546  {
547  private:
548  //==================================== Settings regarding resolutions
549  double mapResolution;
550  unsigned int bandwidth;
551  unsigned int glIntegOrder;
552  unsigned int theta;
553  unsigned int phi;
554  std::vector<unsigned int> bandwidthVec;
555  std::vector<unsigned int> glIntegOrderVec;
556  std::vector<unsigned int> thetaVec;
557  std::vector<unsigned int> phiVec;
558 
559  //==================================== Settings regarding B factors
560  double bFactorValue;
561  double bFactorChange;
562 
563  //==================================== Setting regarding maps and removing noise
564  double noIQRsFromMap;
565 
566  //==================================== Settings regarding concentric shells
567  double shellSpacing;
568  unsigned int manualShells;
569 
570  //==================================== Settings regarding map with phases
571  bool useCOM;
572  bool firstLineCOM;
573 
574  //==================================== Settings regarding space around the structure in lattice
575  double extraSpace;
576  std::vector<double> extraSpaceVec;
577 
578  //==================================== Settings regarding weighting the distances
579  double alpha;
580  double mPower;
581 
582  //==================================== Settings regarding bands to ignore
583  std::vector<int> ignoreLs;
584 
585  //==================================== Settings regarding the structures to use
586  std::vector <std::string> structFiles;
587 
588  //==================================== Settings regarding which distances to compute
589  bool energyLevelDist;
590  bool traceSigmaDist;
591  bool fullRotFnDist;
592 
593  //==================================== Settings regarding dealing with phases
594  bool usePhase;
595  bool saveWithAndWithout;
596 
597  //==================================== Distances holder variables
598  std::vector<double> energyLevelsDistances;
599  std::vector<double> traceSigmaDistances;
600  std::vector<double> fullRotationDistances;
601  std::vector< std::vector<double> > fragEnergyLevelsDistances;
602  std::vector< std::vector<double> > fragTraceSigmaDistances;
603  std::vector< std::vector<double> > fragfullRotationDistances;
604 
605  //==================================== Settings regarding thresholds for hierarchical distance computation
606  double enLevelsThreshold;
607  double trSigmaThreshold;
608 
609  //==================================== Pointer to the comparison object
611 
612  public:
613  //==================================== Public functions
615  ~ProSHADE_distances ( );
616  void saveDatabase ( ProSHADE::ProSHADE_settings* settings );
617  void compareAgainstDatabase ( ProSHADE::ProSHADE_settings* settings,
618  std::vector<std::string>* matchedStrNames );
619  void compareFragAgainstDatabase ( ProSHADE::ProSHADE_settings* settings,
620  std::vector<std::string>* matchedStrNames );
621  std::vector<double> getEnergyLevelsDistances ( );
622  std::vector<double> getTraceSigmaDistances ( );
623  std::vector<double> getFullRotationDistances ( );
624  std::vector< std::vector<double> > getFragEnergyLevelsDistances( );
625  std::vector< std::vector<double> > getFragTraceSigmaDistances ( );
626  std::vector< std::vector<double> > getFragFullRotationDistances( );
627  };
628 
637  {
638  private:
639  //==================================== Private variables
641  double fragBoxSize;
642  std::string mapFragName;
643  double mapFragBoxFraction;
644  double maxMapU;
645  double maxMapV;
646  double maxMapW;
647  double xRange;
648  double yRange;
649  double zRange;
650 
651  double minDensPreNorm;
652  double maxDensPreNorm;
653  double minDensPostNorm;
654  double maxDensPostNorm;
655  double postNormMean;
656  double postNormSdev;
657  double maskVolume;
658  double totalVolume;
659  double maskMean;
660  double maskSdev;
661  double maskMax;
662  double maskMin;
663  double maskDensityRMS;
664  double allDensityRMS;
665  std::array<double,3> origRanges;
666  std::array<double,3> origDims;
667  std::vector<std::string> fragFiles;
668  bool fragFilesExist;
669 
670  public:
671  //==================================== Public functions
674 
675  void printInfo ( int verbose );
676  void fragmentMap ( std::string axOrder, int verbose, ProSHADE::ProSHADE_settings* settings );
677  std::vector<std::string> getFragmentsList ( void );
678  void dealWithHalfMaps ( ProSHADE::ProSHADE_settings* settings, double *xTranslate, double *yTranslate, double *zTranslate );
679  };
680 
681 }
682 
683 
684 //============================================ Include once
685 #endif
This class deals with reading in the data and computing structure specific information including the ...
This class is responsible for reading in map and computing all its features and other possible map ma...
std::vector< std::vector< std::array< double, 6 > > > dnSymmClear
This variable holds the gap corrected Dihedral symmetry results.
unsigned int checkFileType(std::string fileName)
This function checks the input file for being either PDB or MAP formatted.
std::vector< std::array< double, 5 > > cnSymm
This variable holds the complete Cyclic symmetry results.
void printResultsClear(int verbose)
This function prints the cleared results to the screen.
std::vector< std::array< double, 5 > > generateTetrElements(std::vector< std::array< double, 5 > > symmAxes, ProSHADE::ProSHADE_settings *settings, int verbose=0)
This function generates the 12 unique tetrahedral symmetry group elements from the already detected a...
std::vector< std::array< double, 5 > > generateTetrAxes(ProSHADE_comparePairwise *cmpObj, ProSHADE::ProSHADE_settings *settings, std::vector< std::array< double, 5 > > tetrSymm, std::vector< std::array< double, 5 > > allCs, double axisErrorTolerance=0.1, int verbose=0)
This function generates all the tetrahedral symmetry axes from a pair detected by findTetrSymmetry fu...
This namespace contains all the internal objects and their forward declarations.
std::vector< std::array< double, 5 > > icosAxes
This variable holds the 31 unique axes of the octahedral symmetry, if they are found.
std::vector< std::array< double, 5 > > generateOctaElements(std::vector< std::array< double, 5 > > symmAxes, ProSHADE::ProSHADE_settings *settings, int verbose=0)
This function generates the 24 octahedral symmetry group elements from the axes generated by generate...
std::vector< std::array< double, 5 > > generateOctaAxes(ProSHADE_comparePairwise *cmpObj, ProSHADE::ProSHADE_settings *settings, std::vector< std::array< double, 5 > > octaSymm, std::vector< std::array< double, 5 > > allCs, double axisErrorTolerance=0.1, int verbose=0)
This function generates all the ctahedral symmetry elements from a pair detected by findOctaSymmetry ...
This is the executive class responsible for comparing strictly two structures.
std::vector< std::array< double, 5 > > generateIcosAxes(ProSHADE_comparePairwise *cmpObj, ProSHADE::ProSHADE_settings *settings, std::vector< std::array< double, 5 > > icosSymm, std::vector< std::array< double, 5 > > allCs, double axisErrorTolerance=0.1, int verbose=0)
This function generates all the icosahedral symmetry elements from a pair detected by findIcosSymmetr...
std::vector< std::array< double, 5 > > cnSymmClear
This variable holds the gap corrected Cyclic symmetry results.
std::vector< std::array< double, 5 > > generateIcosElements(std::vector< std::array< double, 5 > > symmAxes, ProSHADE::ProSHADE_settings *settings, int verbose=0)
This function generates the 60 icosahedral symmetry group elements from the axes generated by generat...
void printResultsClearHTML(ProSHADE::ProSHADE_settings *settings)
This function prints the cleared results to the HTML file.
bool inputStructureDataType
This variable tells whether input data type is PDB or not.
std::vector< std::array< double, 5 > > getCSymmetries()
This function gives the user programmatical access to the detected C-symmetries.
std::vector< std::array< double, 8 > > getRotFnPeaks()
This function gives the user programmatical access to the symmetry peaks of the symmetry detection pr...
void printResultsRequestHTML(std::string symmetryType, unsigned int symmetryFold, ProSHADE::ProSHADE_settings *settings)
This function prints the cleared results to the HTML report file.
This is the class which computes the symmetry in a single structure.
ProSHADE_symmetry()
Contructor for the ProSHADE_symmetry class.
This is the executive class responsible for comparing two or more structures.
std::vector< std::array< double, 5 > > tetrAxes
This variable holds the 7 unique axes of the tetrahedral symmetry, if they are found.
This is the executive class for computing distances between two or more structures.
This class stores all the settings and is passed to the executive classes instead of multitude of par...
Definition: ProSHADE.h:74
std::vector< std::array< double, 5 > > octaAxes
This variable holds the 13 unique axes of the octahedral symmetry, if they are found.
std::vector< std::vector< std::array< double, 6 > > > getDSymmetries()
This function gives the user programmatical access to the detected D-symmetries.
void printResultsRequest(std::string symmetryType, unsigned int symmetryFold, int verbose)
This function prints the cleared results to the screen.
std::vector< std::vector< std::array< double, 6 > > > dnSymm
This variable holds the complete Dihedral symmetry results.