23 #include <pybind11/pybind11.h>
24 #include <pybind11/stl.h>
25 #include <pybind11/numpy.h>
28 void add_dataClass ( pybind11::module& pyProSHADE )
31 pybind11::class_ < ProSHADE_internal_data::ProSHADE_data > ( pyProSHADE,
"ProSHADE_data" )
34 .def ( pybind11::init ( ) )
35 .def ( pybind11::init ( [] ( std::string strName, pybind11::array_t < float, pybind11::array::c_style | pybind11::array::forcecast > mapData, proshade_single xDmSz, proshade_single yDmSz, proshade_single zDmSz, proshade_unsign xDmInd, proshade_unsign yDmInd, proshade_unsign zDmInd, proshade_signed xFr, proshade_signed yFr, proshade_signed zFr, proshade_signed xT, proshade_signed yT, proshade_signed zT, proshade_unsign inputO )
38 pybind11::buffer_info buf = mapData.request();
39 proshade_unsign len =
static_cast< proshade_unsign
> ( buf.size );
42 double* npVals =
new double[
static_cast<proshade_unsign
> ( len )];
48 for ( proshade_unsign iter = 0; iter < len; iter++ ) { npVals[iter] =
static_cast < double > ( mapData.at(iter) ); }
50 else if ( buf.ndim == 3 )
52 float* dataPtr =
reinterpret_cast < float*
> ( buf.ptr );
53 for ( proshade_unsign xIt = 0; xIt < static_cast<proshade_unsign> ( buf.shape.at(0) ); xIt++ )
55 for ( proshade_unsign yIt = 0; yIt < static_cast<proshade_unsign> ( buf.shape.at(1) ); yIt++ )
57 for ( proshade_unsign zIt = 0; zIt < static_cast<proshade_unsign> ( buf.shape.at(2) ); zIt++ )
59 npVals[zIt +
static_cast< proshade_unsign
> ( buf.shape.at(2) ) * ( yIt +
static_cast< proshade_unsign
> ( buf.shape.at(1) ) * xIt )] =
static_cast < double > ( dataPtr[zIt +
static_cast< proshade_unsign
> ( buf.shape.at(2) ) * ( yIt +
static_cast< proshade_unsign
> ( buf.shape.at(1) ) * xIt )] );
66 std::cerr <<
"!!! ProSHADE PYTHON MODULE ERROR !!! The ProSHADE_data class constructor ( ProSHADE_settings, str, numpy.ndarray, float, float, float, ... ) only supports the third argument input array in the 1D or 3D numpy.ndarray format. The supplied array has " << buf.ndim <<
" dims. Terminating..." << std::endl;
67 exit ( EXIT_FAILURE );
73 static_cast<int> ( len ),
91 .def (
"writeMap", &
ProSHADE_internal_data::ProSHADE_data::writeMap,
"Function for writing out the internal structure representation in MRC MAP format.", pybind11::arg (
"fname" ), pybind11::arg (
"title" ) =
"Created by ProSHADE and written by GEMMI", pybind11::arg (
"mode" ) = 2 )
92 .def (
"writePdb", &
ProSHADE_internal_data::ProSHADE_data::writePdb,
"This function writes out the PDB formatted file coresponding to the structure.", pybind11::arg (
"fname" ), pybind11::arg (
"euA" ) = 0.0, pybind11::arg (
"euB" ) = 0.0, pybind11::arg (
"euG" ) = 0.0, pybind11::arg (
"trsX" ) = 0.0, pybind11::arg (
"trsY" ) = 0.0, pybind11::arg (
"trsZ" ) = 0.0, pybind11::arg (
"firstModel" ) = true )
97 pybind11::array_t < proshade_double > retArr = pybind11::array_t < proshade_double > ( {
self.xDimIndices,
self.yDimIndices,
self.zDimIndices },
98 {
self.yDimIndices *
self.zDimIndices *
sizeof(proshade_double),
99 self.zDimIndices *
sizeof(proshade_double),
100 sizeof(proshade_double) },
108 .def (
"processInternalMap", &
ProSHADE_internal_data::ProSHADE_data::processInternalMap,
"This function simply clusters several map manipulating functions which should be called together. These include centering, phase removal, normalisation, adding extra space, etc.", pybind11::arg (
"settings" ) )
116 .def (
"getReBoxBoundaries",
120 proshade_signed* retVals =
new proshade_signed[6];
126 for ( proshade_unsign iter = 0; iter < 6; iter++ ) { retVals[iter] = settings->
forceBounds[iter]; }
133 static_cast< proshade_signed
> (
self.xDimIndices ),
134 static_cast< proshade_signed
> (
self.yDimIndices ),
135 static_cast< proshade_signed
> (
self.zDimIndices ),
140 self.xDimSize,
self.yDimSize,
self.zDimSize, retVals, settings->
boundsExtraSpace );
146 if ( settings->
useSameBounds && (
self.inputOrder == 0 ) ) {
for ( proshade_unsign iter = 0; iter < 6; iter++ ) { settings->
forceBounds[iter] = retVals[iter]; } }
150 pybind11::capsule pyCapsuleReBox2 ( retVals, [](
void *f ) { proshade_signed* foo =
reinterpret_cast< proshade_signed*
> ( f );
delete foo; } );
153 pybind11::array_t < proshade_signed > retArr = pybind11::array_t < proshade_signed > ( { 6 },
154 {
sizeof(proshade_signed) },
160 },
"This function finds the boundaries enclosing positive map values and adds some extra space." )
161 .def (
"createNewMapFromBounds",
165 proshade_signed* newBounds =
new proshade_signed[6];
169 for ( proshade_unsign iter = 0; iter < 6; iter++ ) { newBounds[iter] = bounds.at(iter); }
172 self.createNewMapFromBounds ( settings, newStr, newBounds );
176 },
"This function creates a new structure from the calling structure and new bounds values." )
193 .def (
"detectSymmetryInStructure",
207 },
"This function runs the symmetry detection algorithms on this structure and saves the results in the settings object.", pybind11::arg (
"settings" ) )
210 .def (
"getRecommendedSymmetryAxes",
214 float* npVals =
new float[
static_cast<unsigned int> ( settings->
detectedSymmetry.size() * 7 )];
218 for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( settings->
detectedSymmetry.size() ); iter++ ) {
for ( proshade_unsign it = 0; it < 7; it++ ) { npVals[(iter*7)+it] =
static_cast< float > ( settings->
detectedSymmetry.at(iter)[it] ); } }
221 pybind11::capsule pyCapsuleStrRecSym ( npVals, [](
void *f ) {
float* foo =
reinterpret_cast< float*
> ( f );
delete foo; } );
224 pybind11::array_t < float > retArr = pybind11::array_t<float> ( {
static_cast<int> ( settings->
detectedSymmetry.size() ),
static_cast<int> ( 7 ) },
225 { 7 *
sizeof(float),
sizeof(
float) },
227 pyCapsuleStrRecSym );
231 },
"This function returns the recommended symmetry axes as a 2D numpy array." )
232 .def (
"getAllCSyms",
236 float* npVals =
new float[
static_cast<unsigned int> ( settings->
allDetectedCAxes.size() * 7 )];
240 for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( settings->
allDetectedCAxes.size() ); iter++ ) {
for ( proshade_unsign it = 0; it < 7; it++ ) { npVals[(iter*7)+it] =
static_cast< float > ( settings->
allDetectedCAxes.at(iter).at(it) ); } }
243 pybind11::capsule pyCapsuleStrSymList ( npVals, [](
void *f ) {
float* foo =
reinterpret_cast< float*
> ( f );
delete foo; } );
246 pybind11::array_t < float > retArr = pybind11::array_t<float> ( {
static_cast<int> ( settings->
allDetectedCAxes.size() ), 7 },
247 { 7 *
sizeof(float),
sizeof(
float) },
249 pyCapsuleStrSymList );
253 },
"This function returns all symmetry axes as a 2D numpy array." )
254 .def (
"getNonCSymmetryAxesIndices",
258 pybind11::dict retDict;
259 pybind11::list dList, tList, oList, iList;
262 for ( proshade_unsign dIt = 0; dIt < static_cast<proshade_unsign> ( settings->
allDetectedDAxes.size() ); dIt++ )
264 pybind11::list memList;
265 for ( proshade_unsign memIt = 0; memIt < static_cast<proshade_unsign> ( settings->
allDetectedDAxes.at(dIt).size() ); memIt++ )
269 dList.append ( memList );
273 for ( proshade_unsign tIt = 0; tIt < static_cast<proshade_unsign> ( settings->
allDetectedTAxes.size() ); tIt++ )
279 for ( proshade_unsign oIt = 0; oIt < static_cast<proshade_unsign> ( settings->
allDetectedOAxes.size() ); oIt++ )
285 for ( proshade_unsign iIt = 0; iIt < static_cast<proshade_unsign> ( settings->
allDetectedIAxes.size() ); iIt++ )
291 retDict[ pybind11::handle ( pybind11::str (
"D" ).ptr ( ) ) ] = dList;
292 retDict[ pybind11::handle ( pybind11::str (
"T" ).ptr ( ) ) ] = tList;
293 retDict[ pybind11::handle ( pybind11::str (
"O" ).ptr ( ) ) ] = oList;
294 retDict[ pybind11::handle ( pybind11::str (
"I" ).ptr ( ) ) ] = iList;
298 },
"This function returns array of non-C axes indices." )
299 .def (
"getAllGroupElements",
303 pybind11::buffer_info buf = axList.request();
304 if ( buf.ndim != 1 ) { std::cerr <<
"!!! ProSHADE PYTHON MODULE ERROR !!! The second argument to getAllGroupElements() must be a 1D numpy array stating the indices of the axes forming the group!" << std::endl; exit ( EXIT_FAILURE ); }
307 std::vector< proshade_unsign > axesList;
308 for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( axList.size() ); iter++ ) {
ProSHADE_internal_misc::addToUnsignVector ( &axesList,
static_cast< proshade_unsign
> ( axList.at(iter) ) ); }
311 std::vector < std::vector< proshade_double > > vals =
self.getAllGroupElements ( settings, axesList, groupType, matrixTolerance );
314 pybind11::list retList;
317 for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( vals.size() ); iter++ )
320 float* npVals =
new float[
static_cast<unsigned int> ( 9 )];
324 for ( proshade_unsign it = 0; it < 9; it++ ) { npVals[it] =
static_cast< float > ( vals.at(iter).at(it) ); }
327 pybind11::capsule pyCapsuleGrpEl ( npVals, [](
void *f ) {
float* foo =
reinterpret_cast< float*
> ( f );
delete foo; } );
330 pybind11::array_t < float > retArr = pybind11::array_t<float> ( { 3, 3 },
331 { 3 *
sizeof(float),
sizeof(
float) },
336 retList.append ( retArr );
341 },
"This function returns the group elements as rotation matrices of any point group described by the detected axes.", pybind11::arg (
"settings" ), pybind11::arg (
"axList" ), pybind11::arg (
"groupType" ) =
"", pybind11::arg(
"matrixTolerance" ) = 0.05 )
343 .def (
"getMapCOMProcessChange",
347 std::vector< proshade_double > vals =
self.getMapCOMProcessChange ();
350 float* npVals =
new float[
static_cast<unsigned int> ( 3 )];
354 for ( proshade_unsign iter = 0; iter < 3; iter++ ) { npVals[iter] =
static_cast< float > ( vals.at(iter) ); }
357 pybind11::capsule pyCapsuleSymShiftDat ( npVals, [](
void *f ) {
float* foo =
reinterpret_cast< float*
> ( f );
delete foo; } );
360 pybind11::array_t < float > retArr = pybind11::array_t<float> ( {
static_cast<int> ( vals.size() ) },
363 pyCapsuleSymShiftDat );
367 },
"This function returns the shift in Angstrom applied to the internal map representation in order to align its COM with the centre of box." )
371 .def (
"getBestRotationMapPeaksEulerAngles",
375 std::vector< proshade_double > vals =
self.getBestRotationMapPeaksEulerAngles ( settings );
378 float* npVals =
new float[
static_cast<unsigned int> (vals.size())];
382 for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( vals.size() ); iter++ ) { npVals[iter] =
static_cast< float > ( vals.at(iter) ); }
385 pybind11::capsule pyCapsuleEuPeak ( npVals, [](
void *f ) {
float* foo =
reinterpret_cast< float*
> ( f );
delete foo; } );
388 pybind11::array_t < float > retArr = pybind11::array_t<float> ( {
static_cast<unsigned int> (vals.size()) },
395 },
"This function returns a vector of three floats, the three Euler angles of the best peak in the rotation map.", pybind11::arg (
"settings" ) )
396 .def (
"getBestRotationMapPeaksRotationMatrix",
400 std::vector< proshade_double > vals =
self.getBestRotationMapPeaksEulerAngles ( settings );
403 proshade_double* retMat =
new proshade_double[9];
408 pybind11::capsule pyCapsuleRMPeak ( retMat, [](
void *f ) { proshade_double* foo =
reinterpret_cast< proshade_double*
> ( f );
delete foo; } );
411 pybind11::array_t < proshade_double > retArr = pybind11::array_t<proshade_double> ( { 3, 3 },
412 { 3 *
sizeof(proshade_double),
sizeof(proshade_double) },
418 },
"This function returns a rotation matrix representing the best peak in the rotation map.", pybind11::arg (
"settings" ) )
423 .def (
"getOverlayTranslations",
427 std::vector< proshade_double > vals =
self.getBestTranslationMapPeaksAngstrom ( staticStructure, eulA, eulB, eulG );
430 pybind11::dict retDict;
431 pybind11::list rotCen, toOverlay;
434 rotCen.append (
self.originalPdbRotCenX );
435 rotCen.append (
self.originalPdbRotCenY );
436 rotCen.append (
self.originalPdbRotCenZ );
438 toOverlay.append (
self.originalPdbTransX );
439 toOverlay.append (
self.originalPdbTransY );
440 toOverlay.append (
self.originalPdbTransZ );
443 retDict[ pybind11::handle ( pybind11::str (
"centreOfRotation" ).ptr ( ) ) ] = rotCen;
444 retDict[ pybind11::handle ( pybind11::str (
"rotCenToOverlay" ).ptr ( ) ) ] = toOverlay;
448 },
"This function returns the vector from optimal rotation centre to origin and the optimal overlay translation vector. These two vectors allow overlaying the inputs (see documentation for details on how the two vectors should be used).", pybind11::arg (
"staticStructure" ), pybind11::arg (
"eulA" ), pybind11::arg (
"eulB" ), pybind11::arg (
"eulG" ) )
449 .def (
"translateMap", &
ProSHADE_internal_data::ProSHADE_data::translateMap,
"This function translates the map by a given number of Angstroms along the three axes. Please note the translation happens firstly to the whole map box and only the translation remainder that cannot be achieved by moving the box will be corrected for using reciprocal space translation within the box.", pybind11::arg (
"trsX" ), pybind11::arg (
"trsY" ), pybind11::arg (
"trsZ" ) )
452 .def (
"findSHIndex",
456 proshade_signed index = seanindex (
static_cast< int > ( order ),
457 static_cast< int > ( band ),
458 static_cast< int > (
self.spheres[shell]->getLocalBandwidth() ) );
462 },
"This function finds the correct index for given shell, band and order in the spherical harmonics array. Please note that the order is expected in range -band <= 0 <- band and NOT from 0 to ( 2 * band ) + 1." )
463 .def (
"getSphericalHarmonics",
467 std::complex<proshade_double>* npVals =
new std::complex<proshade_double>[
static_cast<proshade_unsign
> (
self.noSpheres * pow(
self.maxShellBand, 2.0 ) )];
471 for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> (
self.noSpheres * pow(
self.maxShellBand, 2.0 ) ); iter++ ) { npVals[iter].real ( 0.0 ); npVals[iter].imag ( 0.0 ); }
474 proshade_signed pyPosSH;
475 proshade_signed pyPos;
478 for ( proshade_signed shIt = 0; shIt < static_cast<proshade_signed> (
self.noSpheres ); shIt++ )
480 for ( proshade_signed bnd = 0; bnd < static_cast<proshade_signed> (
self.spheres[shIt]->getLocalBandwidth() ); bnd++ )
482 for ( proshade_signed order = -bnd; order <= bnd; order++ )
484 pyPosSH = (
static_cast< proshade_signed
> ( shIt ) *
static_cast< proshade_signed
> ( std::pow (
self.maxShellBand, 2 ) ) );
485 pyPos = seanindex (
static_cast< int > ( order ),
486 static_cast< int > ( bnd ),
487 static_cast< int > (
self.spheres[shIt]->getLocalBandwidth() ) );
488 npVals[pyPosSH+pyPos].real (
self.sphericalHarmonics[shIt][pyPos][0] );
489 npVals[pyPosSH+pyPos].imag (
self.sphericalHarmonics[shIt][pyPos][1] );
495 pybind11::capsule pyCapsuleSHs ( npVals, [](
void *f ) { std::complex<proshade_double>* foo =
reinterpret_cast< std::complex<proshade_double>*
> ( f );
delete foo; } );
498 pybind11::array_t < std::complex<proshade_double> > retArr = pybind11::array_t < std::complex<proshade_double > > (
499 {
static_cast<int> (
self.noSpheres ),
static_cast<int> ( pow (
self.maxShellBand, 2.0 ) ) },
500 {
sizeof ( std::complex < proshade_double > ) *
static_cast< proshade_unsign
> ( std::pow (
static_cast< proshade_unsign
> (
self.maxShellBand ), 2 ) ),
sizeof ( std::complex < proshade_double > ) },
506 },
"This function returns a 2D numpy array of complex numbers representing the spherical harmonics computed for the structure. The first dimension of the array is the spheres (i.e. each sphere has its own array) and the second dimension is the band and order combination as given by the findSHIndex() function. Please note that while each sphere can have different number of spherical harmonics coefficients (depending on the settings.progressiveSphereMapping value), this array uses maxShellBand**2 values to make sure the length are equal for all spheres. To avoid issues, please use the findSHIndex() to correctly find the index of a particular shell, band, order combination." )
511 std::complex<proshade_double>* npVals =
new std::complex < proshade_double >[
static_cast<proshade_unsign
> (
self.maxShellBand * ( (
self.maxShellBand * 2 ) + 1 ) * ( (
self.maxShellBand * 2 ) + 1 ) )];
515 proshade_signed index = 0;
518 for ( proshade_unsign iter = 0; iter < static_cast < proshade_unsign > (
self.maxShellBand * ( (
self.maxShellBand * 2 ) + 1 ) * ( (
self.maxShellBand * 2 ) + 1 ) ); iter++ ) { npVals[iter].real ( 0.0 ); npVals[iter].imag ( 0.0 ); }
521 for ( proshade_signed bandIter = 0; bandIter < static_cast< proshade_signed > (
self.maxShellBand ); bandIter++ )
523 for ( proshade_signed order1 = 0; order1 < ( ( bandIter * 2 ) + 1 ); order1++ )
525 for ( proshade_signed order2 = 0; order2 < ( ( bandIter * 2 ) + 1 ); order2++ )
527 index = order2 + ( (
static_cast< proshade_signed
> (
self.maxShellBand ) * 2 ) + 1 ) * ( order1 + ( (
static_cast< proshade_signed
> (
self.maxShellBand ) * 2 ) + 1 ) * bandIter );
528 npVals[index].real (
self.eMatrices[bandIter][order1][order2][0] );
529 npVals[index].imag (
self.eMatrices[bandIter][order1][order2][1] );
535 pybind11::capsule pyCapsuleEMs ( npVals, [](
void *f ) { std::complex<proshade_double>* foo =
reinterpret_cast< std::complex<proshade_double>*
> ( f );
delete foo; } );
538 pybind11::array_t < std::complex<proshade_double> > retArr = pybind11::array_t < std::complex<proshade_double > > (
539 {
self.maxShellBand, ( (
self.maxShellBand * 2 ) + 1 ), ( (
self.maxShellBand * 2 ) + 1 ) },
540 {
sizeof ( std::complex < proshade_double > ) *
static_cast< proshade_unsign
> ( ( (
static_cast< proshade_unsign
> (
self.maxShellBand ) * 2 ) + 1 ) * ( (
static_cast< proshade_unsign
> (
self.maxShellBand ) * 2 ) + 1 ) ),
541 sizeof ( std::complex < proshade_double > ) * (
static_cast< proshade_unsign
> (
self.maxShellBand ) * 2 ) + 1,
542 sizeof ( std::complex < proshade_double > ) },
548 },
"This function returns the E matrix values (these are the integral over all spheres of ( c1^(l,m) * c2^(l,m') ) values) obtained when rotation function or self-rotation function are computed. The returned array has three dimensions, first being the band, second being the order1 and third being the order2. Please note that as negative indexing does not simply work, the order indexing starts from 0 - i.e. array[1][0][0] means band 1 ; order1 = -1 and order2 = -1." )
549 .def (
"getSO3Coefficients",
553 std::complex<proshade_double>* npVals =
new std::complex < proshade_double >[
static_cast<proshade_unsign
> (
self.maxShellBand * ( (
self.maxShellBand * 2 ) + 1 ) * ( (
self.maxShellBand * 2 ) + 1 ) )];
557 proshade_signed index = 0;
560 for ( proshade_unsign iter = 0; iter < static_cast < proshade_unsign > (
self.maxShellBand * ( (
self.maxShellBand * 2 ) + 1 ) * ( (
self.maxShellBand * 2 ) + 1 ) ); iter++ ) { npVals[iter].real ( 0.0 ); npVals[iter].imag ( 0.0 ); }
563 for ( proshade_signed bandIter = 0; bandIter < static_cast< proshade_signed > (
self.maxShellBand ); bandIter++ )
565 for ( proshade_signed order1 = 0; order1 < ( ( bandIter * 2 ) + 1 ); order1++ )
567 for ( proshade_signed order2 = 0; order2 < ( ( bandIter * 2 ) + 1 ); order2++ )
569 index = order2 + ( (
static_cast< proshade_signed
> (
self.maxShellBand ) * 2 ) + 1 ) * ( order1 + ( (
static_cast< proshade_signed
> (
self.maxShellBand ) * 2 ) + 1 ) * bandIter );
570 npVals[index].real (
self.so3Coeffs[
self.so3CoeffsArrayIndex (
static_cast< proshade_signed
> ( order1 - bandIter ), order2 - bandIter, bandIter )][0] );
571 npVals[index].imag (
self.so3Coeffs[
self.so3CoeffsArrayIndex (
static_cast< proshade_signed
> ( order1 - bandIter ), order2 - bandIter, bandIter )][1] );
577 pybind11::capsule pyCapsuleSOCoeffs ( npVals, [](
void *f ) { std::complex<proshade_double>* foo =
reinterpret_cast< std::complex<proshade_double>*
> ( f );
delete foo; } );
580 pybind11::array_t < std::complex<proshade_double> > retArr = pybind11::array_t < std::complex<proshade_double > > (
581 {
self.maxShellBand, ( (
self.maxShellBand * 2 ) + 1 ), ( (
self.maxShellBand * 2 ) + 1 ) },
582 {
sizeof ( std::complex < proshade_double > ) * ( (
static_cast< proshade_unsign
> (
self.maxShellBand ) * 2 ) + 1 ) * ( (
static_cast< proshade_unsign
> (
self.maxShellBand ) * 2 ) + 1 ),
583 sizeof ( std::complex < proshade_double > ) * (
static_cast< proshade_unsign
> (
self.maxShellBand ) * 2 ) + 1,
584 sizeof ( std::complex < proshade_double > ) },
590 },
"This function returns the SO(3) coefficient values (these are normalised E matrix values with sign modification) obtained when rotation function or self-rotation function are computed. The returned array has three dimensions, first being the band, second being the order1 and third being the order2. Please note that as negative indexing does not simply work, the order indexing starts from 0 - i.e. array[1][0][0] means band 1 ; order1 = -1 and order2 = -1." )
591 .def (
"getRotationFunctionMap",
595 std::complex<proshade_double>* npVals =
new std::complex < proshade_double >[
static_cast<proshade_unsign
> ( (
self.maxShellBand * 2 ) * (
self.maxShellBand * 2 ) * (
self.maxShellBand * 2 ) )];
599 for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( (
self.maxShellBand * 2 ) * (
self.maxShellBand * 2 ) * (
self.maxShellBand * 2 ) ); iter++ ) { npVals[iter].real (
self.so3CoeffsInverse[iter][0] ); npVals[iter].imag (
self.so3CoeffsInverse[iter][1] ); }
602 pybind11::capsule pyCapsuleRotMap ( npVals, [](
void *f ) { std::complex<proshade_double>* foo =
reinterpret_cast< std::complex<proshade_double>*
> ( f );
delete foo; } );
605 pybind11::array_t < std::complex < proshade_double > > retArr = pybind11::array_t < std::complex < proshade_double > > (
606 { (
self.maxShellBand * 2 ), (
self.maxShellBand * 2 ), (
self.maxShellBand * 2 ) },
607 { (
self.maxShellBand * 2 ) * (
self.maxShellBand * 2 ) *
sizeof(std::complex < proshade_double >),
608 (
self.maxShellBand * 2 ) *
sizeof(std::complex < proshade_double >),
609 sizeof(std::complex < proshade_double >) },
615 },
"This function returns the (self) rotation function as a three-dimensional map of complex numbers." )
616 .def (
"getRotationMatrixFromSOFTCoordinates",
620 proshade_double* npVals =
new proshade_double[9];
624 proshade_double eulA, eulB, eulG;
633 pybind11::capsule pyCapsuleRMSoft ( npVals, [](
void *f ) { proshade_double* foo =
reinterpret_cast< proshade_double*
> ( f );
delete foo; } );
636 pybind11::array_t < proshade_double > retArr = pybind11::array_t < proshade_double > ( { 3, 3 },
637 { 3 *
sizeof(proshade_double),
sizeof(proshade_double) },
643 },
"This function converts a given rotation function map position onto the corresponding rotation matrix.", pybind11::arg (
"xPos" ), pybind11::arg (
"yPos" ), pybind11::arg (
"zPos" ) )
644 .def (
"getTranslationFunctionMap",
648 std::complex<proshade_double>* npVals =
new std::complex < proshade_double >[
static_cast<proshade_unsign
> (
self.getXDim() *
self.getYDim() *
self.getZDim() )];
652 for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> (
self.getXDim() *
self.getYDim() *
self.getZDim() ); iter++ ) { npVals[iter].real (
self.translationMap[iter][0] ); npVals[iter].imag (
self.translationMap[iter][1] ); }
655 pybind11::capsule pyCapsuleTrsMap ( npVals, [](
void *f ) { std::complex<proshade_double>* foo =
reinterpret_cast< std::complex<proshade_double>*
> ( f );
delete foo; } );
658 pybind11::array_t < std::complex < proshade_double > > retArr = pybind11::array_t < std::complex < proshade_double > > (
659 {
self.getXDim(),
self.getYDim(),
self.getZDim() },
660 {
self.getYDim() *
self.getZDim() *
sizeof(std::complex < proshade_double >),
661 self.getZDim() *
sizeof(std::complex < proshade_double >),
662 sizeof(std::complex < proshade_double >) },
668 },
"This function returns the translation function as a three-dimensional map of complex numbers." )
706 .def (
"__repr__", [] ( ) {
return "<ProSHADE_data class object> (This class contains all information, results and available functionalities for a structure)"; } );
709 pyProSHADE.def (
"computeGroupElementsForGroup",
710 [] ( proshade_double x, proshade_double y, proshade_double z, proshade_unsign fold ) -> pybind11::array_t < proshade_double >
716 proshade_double* npVals =
new proshade_double[
static_cast<proshade_unsign
> ( retVec.size() ) * 9];
720 for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( retVec.size() ); iter++ ) {
for ( proshade_unsign it = 0; it < 9; it++ ) { npVals[(iter*9)+it] = retVec.at(iter).at(it); } }
723 pybind11::capsule pyCapsuleGrEls ( npVals, [](
void *f ) { proshade_double* foo =
reinterpret_cast< proshade_double*
> ( f );
delete foo; } );
726 pybind11::array_t < proshade_double > retArr = pybind11::array_t < proshade_double > (
727 {
static_cast< int > ( retVec.size() ), 3, 3 },
728 { 9 *
sizeof(proshade_double), 3 *
sizeof(proshade_double),
sizeof(proshade_double) },
734 },
"This function takes the axis and fold and computes all resulting group elements (as rotation matrices), returning them in numpy.ndarray.", pybind11::arg (
"x" ), pybind11::arg (
"y" ), pybind11::arg (
"z" ), pybind11::arg (
"fold" ) );
735 pyProSHADE.def (
"joinElementsFromDifferentGroups",
736 [] ( pybind11::array_t < proshade_double, pybind11::array::c_style | pybind11::array::forcecast > first,
737 pybind11::array_t < proshade_double, pybind11::array::c_style | pybind11::array::forcecast > second,
738 proshade_double matrixTolerance,
739 bool combine ) -> pybind11::array_t < proshade_double >
742 pybind11::buffer_info buf1 = first.request();
743 pybind11::buffer_info buf2 = second.request();
746 if ( buf1.ndim != 3 || buf2.ndim != 3 )
748 std::cerr <<
"Function joinElementsFromDifferentGroups() arguments first and second should be numpy.ndarrays with 3 dimensions indexed as follos: first[elementNumber][elementRotationMatrixRow][elementRotationMatrixColumn] - the same format as returned by the computeGroupElementsForGroup() function." << std::endl;
749 return ( pybind11::array_t < proshade_double > () );
753 std::vector< std::vector< proshade_double > > fVec, sVec;
755 proshade_double* dataPtr1 =
reinterpret_cast < proshade_double*
> ( buf1.ptr );
756 for ( proshade_unsign elIt = 0; elIt < static_cast<proshade_unsign> ( buf1.shape.at(0) ); elIt++ )
758 std::vector< proshade_double > rotMat;
759 for ( proshade_unsign rowIt = 0; rowIt < static_cast<proshade_unsign> ( buf1.shape.at(1) ); rowIt++ )
761 for ( proshade_unsign colIt = 0; colIt < static_cast<proshade_unsign> ( buf1.shape.at(2) ); colIt++ )
769 proshade_double* dataPtr2 =
reinterpret_cast < proshade_double*
> ( buf2.ptr );
770 for ( proshade_unsign elIt = 0; elIt < static_cast<proshade_unsign> ( buf2.shape.at(0) ); elIt++ )
772 std::vector< proshade_double > rotMat;
773 for ( proshade_unsign rowIt = 0; rowIt < static_cast<proshade_unsign> ( buf2.shape.at(1) ); rowIt++ )
775 for ( proshade_unsign colIt = 0; colIt < static_cast<proshade_unsign> ( buf2.shape.at(2) ); colIt++ )
788 proshade_double* npVals =
new proshade_double[
static_cast<proshade_unsign
> ( retVec.size() ) * 9];
792 for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( retVec.size() ); iter++ ) {
for ( proshade_unsign it = 0; it < 9; it++ ) { npVals[(iter*9)+it] = retVec.at(iter).at(it); } }
795 pybind11::capsule pyCapsuleGrElsCombo ( npVals, [](
void *f ) { proshade_double* foo =
reinterpret_cast< proshade_double*
> ( f );
delete foo; } );
798 pybind11::array_t < proshade_double > retArr = pybind11::array_t < proshade_double > (
799 {
static_cast< int > ( retVec.size() ), 3, 3 },
800 { 9 *
sizeof(proshade_double), 3 *
sizeof(proshade_double),
sizeof(proshade_double) },
802 pyCapsuleGrElsCombo );
806 },
"This function takes the axis and fold and computes all resulting group elements (as rotation matrices), returning them in numpy.ndarray.", pybind11::arg (
"first" ), pybind11::arg (
"second" ), pybind11::arg (
"matrixTolerance" ), pybind11::arg (
"combine" ) );