ProSHADE  0.7.6.6 (JUL 2022)
Protein Shape Detection
pyProSHADE_symmetry.cpp
Go to the documentation of this file.
1 
22 //==================================================== Include PyBind11 header
23 #include <pybind11/pybind11.h>
24 #include <pybind11/stl.h>
25 #include <pybind11/numpy.h>
26 
27 //==================================================== Add the ProSHADE_settings and ProSHADE_run classes to the PyBind11 module
28 void add_symmetryNamespace ( pybind11::module& pyProSHADE )
29 {
30 // pyProSHADE.def ( "findReliableUnphasedSymmetries",
31 // [] ( ProSHADE_settings* settings, proshade_signed verbose, proshade_signed messageShift, proshade_double tolerance ) -> pybind11::array_t < proshade_unsign >
32 // {
33 // //== Run the detection
34 // std::vector< proshade_unsign > rels = ProSHADE_internal_symmetry::findReliableUnphasedSymmetries ( &settings->allDetectedCAxes, verbose, messageShift, tolerance );
35 //
36 // //== Allocate memory for the numpy values
37 // proshade_unsign* npVals = new proshade_unsign[static_cast<unsigned int> ( rels.size() )];
38 // ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
39 //
40 // //== Copy values
41 // for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( rels.size() ); iter++ ) { npVals[iter] = static_cast< proshade_unsign > ( rels.at(iter) ); }
42 //
43 // //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
44 // pybind11::capsule pyCapsuleRelSyms ( npVals, []( void *f ) { proshade_unsign* foo = reinterpret_cast< proshade_unsign* > ( f ); delete foo; } );
45 //
46 // //== Copy the value
47 // pybind11::array_t < proshade_unsign > retArr = pybind11::array_t< proshade_unsign > ( { static_cast<int> ( rels.size() ) }, // Shape
48 // { sizeof(proshade_unsign) }, // C-stype strides
49 // npVals, // Data
50 // pyCapsuleRelSyms ); // Capsule
51 //
52 // //== Done
53 // return ( retArr );
54 // }, "This function checks the list of detected axes (presumably from phaseless symmetry detection) and returns the best dihedral (or cyclic, if no dihedral is found) point group, or empty vector if nothing is found.", pybind11::arg ( "settings" ), pybind11::arg ( "verbose" ), pybind11::arg ( "messageShift" ) = 1, pybind11::arg ( "tolerance" ) = 0.1 );
55 
56  pyProSHADE.def ( "optimiseDGroupAngleFromAxesHeights",
57  [] ( pybind11::array_t < proshade_unsign > selection, ProSHADE_internal_data::ProSHADE_data* dataObj, ProSHADE_settings* settings )
58  {
59  //== Sanity check
60  pybind11::buffer_info selection_buf = selection.request();
61  if ( selection_buf.ndim != 1 ) { std::cerr << "!!! ProSHADE PYTHON MODULE ERROR !!! The first argument to optimiseDGroupAngleFromAxesHeights() must be a 1D numpy array!" << std::endl; exit ( EXIT_FAILURE ); }
62  if ( selection_buf.shape.at(0) != 2 ) { std::cerr << "!!! ProSHADE PYTHON MODULE ERROR !!! The first argument to optimiseDGroupAngleFromAxesHeights() must be array of length 2!" << std::endl; exit ( EXIT_FAILURE ); }
63 
64  //== Copy data to C++ format
65  proshade_unsign* arrStart = static_cast< proshade_unsign* > ( selection_buf.ptr );
66  std::vector< proshade_unsign > sel;
67  ProSHADE_internal_misc::addToUnsignVector ( &sel, arrStart[0] );
68  ProSHADE_internal_misc::addToUnsignVector ( &sel, arrStart[1] );
69 
70  //== Convert data type
71  std::vector< std::vector< proshade_double > > allCs;
72  std::vector< proshade_double > hlpVec;
73  for ( size_t it1 = 0; it1 < dataObj->getCyclicAxes()->size(); it1++ )
74  {
75  hlpVec.clear ( );
76  for ( size_t it2 = 0; it2 < 7; it2++ ) { ProSHADE_internal_misc::addToDoubleVector ( &hlpVec, dataObj->getCyclicAxes()->at(it1)[it2] ); }
78  }
79 
80  //== Call the C++ function
81  ProSHADE_internal_symmetry::optimiseDGroupAngleFromAxesHeights ( &allCs, sel, dataObj, settings );
82 
83  //== Done (all changes are in the data object, so nothing needs to be passed to Python)
84  return ;
85  }, "This function takes two axes with almost dihedral angle and optimises their relative positions as well as orientation with respect to the optimal angle and the rotation function.", pybind11::arg ( "selection" ), pybind11::arg ( "dataObj" ), pybind11::arg ( "settings" ) );
86 
87  pyProSHADE.def ( "findPointFromTranslations",
88  [] ( ProSHADE_settings* settings, ProSHADE_internal_data::ProSHADE_data* dataObj, pybind11::array_t < float > allCs, proshade_unsign whichAxis ) -> pybind11::array_t < proshade_double >
89  {
90  //== Sanity check
91  pybind11::buffer_info allCs_buf = allCs.request();
92  if ( allCs_buf.ndim != 2 ) { std::cerr << "!!! ProSHADE PYTHON MODULE ERROR !!! The third argument to findPointFromTranslations() must be a 2D numpy array!" << std::endl; exit ( EXIT_FAILURE ); }
93 
94  //== Copy to C++ data format
95  std::vector< std::vector < proshade_double > > allCAxes;
96  std::vector < proshade_double > hlpVec;
97  float* arrStart = static_cast< float* > ( allCs_buf.ptr );
98  for ( proshade_unsign axIt = 0; axIt < static_cast< proshade_unsign > ( allCs_buf.shape.at(0) ); axIt++ )
99  {
100  //== Clean up
101  hlpVec.clear ( );
102 
103  //== Copy vals
104  for ( proshade_unsign vIt = 0; vIt < 7; vIt++ )
105  {
106  ProSHADE_internal_misc::addToDoubleVector ( &hlpVec, static_cast< proshade_double > ( arrStart[(axIt*7)+vIt] ) );
107  }
108 
110  }
111 
112  //== Prepare all required memory
113  fftw_complex *origMap = nullptr, *origCoeffs = nullptr, *rotMapComplex = nullptr, *rotCoeffs = nullptr, *trFunc = nullptr, *trFuncCoeffs = nullptr;
114  fftw_plan planForwardFourier, planForwardFourierRot, planReverseFourierComb;
115  ProSHADE_internal_symmetry::allocateCentreOfMapFourierTransforms ( dataObj->getXDim(), dataObj->getYDim(), dataObj->getZDim(), origMap, origCoeffs, rotMapComplex, rotCoeffs, trFunc, trFuncCoeffs, &planForwardFourier, &planForwardFourierRot, &planReverseFourierComb );
116 
117  //== Compute Fourier for the original map
118  for ( size_t it = 0; it < static_cast< size_t > ( dataObj->getXDim() * dataObj->getYDim() * dataObj->getZDim() ); it++ ) { origMap[it][0] = dataObj->getMapValue( it ); origMap[it][1] = 0.0; }
119  fftw_execute ( planForwardFourier );
120 
121  //== Run C++ code
122  std::vector< proshade_unsign > axLst;
123  std::vector< std::vector< proshade_double > > symElems;
124  ProSHADE_internal_misc::addToUnsignVector ( &axLst, whichAxis );
125  symElems = dataObj->getAllGroupElements ( &allCAxes, axLst, "C", settings->axisErrTolerance );
126  std::vector< proshade_double > pointPos = ProSHADE_internal_symmetry::findPointFromTranslations ( dataObj,
127  symElems,
128  origCoeffs, rotMapComplex,
129  rotCoeffs, planForwardFourierRot,
130  trFuncCoeffs, trFunc,
131  planReverseFourierComb );
132 
133  //== Release required memory
134  ProSHADE_internal_symmetry::releaseCentreOfMapFourierTransforms ( origMap, origCoeffs, rotMapComplex, rotCoeffs, trFunc, trFuncCoeffs, planForwardFourier, planForwardFourierRot, planReverseFourierComb );
135 
136  //== Allocate memory for the numpy values
137  proshade_double* npVals = new proshade_double[3];
138  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
139 
140  //== Copy values
141  for ( proshade_unsign iter = 0; iter < 3; iter++ ) { npVals[iter] = static_cast< proshade_double > ( pointPos.at(iter) ); }
142 
143  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
144  pybind11::capsule pyCapsuleSymPoint ( npVals, []( void *f ) { proshade_double* foo = reinterpret_cast< proshade_double* > ( f ); delete foo; } );
145 
146  //== Copy the value
147  pybind11::array_t < proshade_double > retArr = pybind11::array_t< proshade_double > ( { static_cast<int> ( 3 ) }, // Shape
148  { sizeof(proshade_double) }, // C-stype strides
149  npVals, // Data
150  pyCapsuleSymPoint ); // Capsule
151 
152  //== Done
153  return ( retArr );
154 
155 
156  }, "This function computes the average of optimal translations for a cyclic point group.", pybind11::arg ( "settings" ), pybind11::arg ( "dataObj" ), pybind11::arg ( "allCs" ), pybind11::arg ( "whichAxis" ) );
157 }
ProSHADE_internal_symmetry::findPointFromTranslations
std::vector< proshade_double > findPointFromTranslations(ProSHADE_internal_data::ProSHADE_data *symStr, std::vector< std::vector< proshade_double > > symElems, fftw_complex *origCoeffs, fftw_complex *rotMapComplex, fftw_complex *rotCoeffs, fftw_plan planForwardFourierRot, fftw_complex *trFuncCoeffs, fftw_complex *trFunc, fftw_plan planReverseFourierComb)
This function computes the average of optimal translations for a cyclic point group.
Definition: ProSHADE_symmetry.cpp:3985
ProSHADE_internal_symmetry::allocateCentreOfMapFourierTransforms
void allocateCentreOfMapFourierTransforms(proshade_unsign xDim, proshade_unsign yDim, proshade_unsign zDim, fftw_complex *&origMap, fftw_complex *&origCoeffs, fftw_complex *&rotMapComplex, fftw_complex *&rotCoeffs, fftw_complex *&trFunc, fftw_complex *&trFuncCoeffs, fftw_plan *planForwardFourier, fftw_plan *planForwardFourierRot, fftw_plan *planReverseFourierComb)
This function allocates the required memory for the Fourier transforms required to find the centre of...
Definition: ProSHADE_symmetry.cpp:3841
ProSHADE_internal_data::ProSHADE_data
This class contains all inputed and derived data for a single structure.
Definition: ProSHADE_data.hpp:49
ProSHADE_internal_symmetry::optimiseDGroupAngleFromAxesHeights
void optimiseDGroupAngleFromAxesHeights(std::vector< std::vector< proshade_double > > *ret, ProSHADE_internal_data::ProSHADE_data *dataObj, ProSHADE_settings *settings)
This function takes two axes with almost dihedral angle and optimises their relative positions as wel...
Definition: ProSHADE_symmetry.cpp:3269
ProSHADE_internal_data::ProSHADE_data::getAllGroupElements
std::vector< std::vector< proshade_double > > getAllGroupElements(std::vector< proshade_unsign > axesList, std::string groupType="", proshade_double matrixTolerance=0.05)
This function returns the group elements as rotation matrices of any defined point group.
Definition: ProSHADE_data.cpp:3122
ProSHADE_internal_misc::addToDoubleVector
void addToDoubleVector(std::vector< proshade_double > *vecToAddTo, proshade_double elementToAdd)
Adds the element to the vector.
Definition: ProSHADE_misc.cpp:77
ProSHADE_internal_data::ProSHADE_data::getXDim
proshade_unsign getXDim(void)
This function allows access to the map size in indices along the X axis.
Definition: ProSHADE_data.cpp:4074
ProSHADE_internal_data::ProSHADE_data::getCyclicAxes
std::vector< proshade_double * > * getCyclicAxes(void)
This function allows access to the list of detected cyclic axes.
Definition: ProSHADE_data.cpp:4239
ProSHADE_internal_data::ProSHADE_data::getYDim
proshade_unsign getYDim(void)
This function allows access to the map size in indices along the Y axis.
Definition: ProSHADE_data.cpp:4084
ProSHADE_settings
This class stores all the settings and is passed to the executive classes instead of a multitude of p...
Definition: ProSHADE_settings.hpp:37
ProSHADE_internal_data::ProSHADE_data::getZDim
proshade_unsign getZDim(void)
This function allows access to the map size in indices along the Z axis.
Definition: ProSHADE_data.cpp:4094
ProSHADE_internal_data::ProSHADE_data::getMapValue
proshade_double getMapValue(proshade_unsign pos)
This function returns the internal map representation value of a particular array position.
Definition: ProSHADE_data.cpp:3738
ProSHADE_internal_symmetry::releaseCentreOfMapFourierTransforms
void releaseCentreOfMapFourierTransforms(fftw_complex *origMap, fftw_complex *origCoeffs, fftw_complex *rotMapComplex, fftw_complex *rotCoeffs, fftw_complex *trFunc, fftw_complex *trFuncCoeffs, fftw_plan planForwardFourier, fftw_plan planForwardFourierRot, fftw_plan planReverseFourierComb)
This function releases the allocated memory for the Fourier transforms used to find the centre of the...
Definition: ProSHADE_symmetry.cpp:3881
ProSHADE_internal_misc::checkMemoryAllocation
void checkMemoryAllocation(chVar checkVar, std::string fileP, unsigned int lineP, std::string funcP, std::string infoP="This error may occurs when ProSHADE requests memory to be\n : allocated to it and this operation fails. This could\n : happen when not enough memory is available, either due to\n : other processes using a lot of memory, or when the machine\n : does not have sufficient memory available. Re-run to see\n : if this problem persists.")
Checks if memory was allocated properly.
Definition: ProSHADE_misc.hpp:73
ProSHADE_internal_misc::addToDoubleVectorVector
void addToDoubleVectorVector(std::vector< std::vector< proshade_double > > *vecToAddTo, std::vector< proshade_double > elementToAdd)
Adds the element to the vector of vectors.
Definition: ProSHADE_misc.cpp:233
ProSHADE_internal_misc::addToUnsignVector
void addToUnsignVector(std::vector< proshade_unsign > *vecToAddTo, proshade_unsign elementToAdd)
Adds the element to the vector.
Definition: ProSHADE_misc.cpp:99
ProSHADE_settings::axisErrTolerance
proshade_double axisErrTolerance
Allowed error on vector axis in in dot product ( acos ( 1 - axErr ) is the allowed difference in radi...
Definition: ProSHADE_settings.hpp:132