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 < ProSHADE_settings* > ( ), pybind11::arg (
"settings" ) )
35 .def ( pybind11::init ( [] (
ProSHADE_settings* settings, 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 = 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 + buf.shape.at(2) * ( yIt + buf.shape.at(1) * xIt )] =
static_cast < double > ( dataPtr[zIt + buf.shape.at(2) * ( yIt + 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 );
74 static_cast<int> ( len ),
92 .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 )
93 .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 )
98 pybind11::array_t < proshade_double > retArr = pybind11::array_t < proshade_double > ( {
self.xDimIndices,
self.yDimIndices,
self.zDimIndices },
99 {
self.yDimIndices *
self.zDimIndices *
sizeof(proshade_double),
100 self.zDimIndices *
sizeof(proshade_double),
101 sizeof(proshade_double) },
109 .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" ) )
117 .def (
"getReBoxBoundaries",
121 proshade_signed* retVals =
new proshade_signed[6];
127 for ( proshade_unsign iter = 0; iter < 6; iter++ ) { retVals[iter] = settings->
forceBounds[iter]; }
134 self.xDimSize,
self.yDimSize,
self.zDimSize, retVals );
138 self.xDimSize,
self.yDimSize,
self.zDimSize, retVals, settings->
boundsExtraSpace );
144 if ( settings->
useSameBounds && (
self.inputOrder == 0 ) ) {
for ( proshade_unsign iter = 0; iter < 6; iter++ ) { settings->
forceBounds[iter] = retVals[iter]; } }
148 pybind11::capsule pyCapsuleReBox2 ( retVals, [](
void *f ) { proshade_signed* foo =
reinterpret_cast< proshade_signed*
> ( f );
delete foo; } );
151 pybind11::array_t < proshade_signed > retArr = pybind11::array_t < proshade_signed > ( { 6 },
152 {
sizeof(proshade_signed) },
158 },
"This function finds the boundaries enclosing positive map values and adds some extra space." )
159 .def (
"createNewMapFromBounds",
163 proshade_signed* newBounds =
new proshade_signed[6];
167 for ( proshade_unsign iter = 0; iter < 6; iter++ ) { newBounds[iter] = bounds.at(iter); }
170 self.createNewMapFromBounds ( settings, newStr, newBounds );
174 },
"This function creates a new structure from the calling structure and new bounds values." )
191 .def (
"detectSymmetryInStructure",
205 },
"This function runs the symmetry detection algorithms on this structure and saves the results in the settings object.", pybind11::arg (
"settings" ) )
208 .def (
"getRecommendedSymmetryAxes",
212 float* npVals =
new float[
static_cast<unsigned int> ( settings->
detectedSymmetry.size() * 6 )];
216 for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( settings->
detectedSymmetry.size() ); iter++ ) {
for ( proshade_unsign it = 0; it < 6; it++ ) { npVals[(iter*6)+it] = settings->
detectedSymmetry.at(iter)[it]; } }
219 pybind11::capsule pyCapsuleStrRecSym ( npVals, [](
void *f ) {
float* foo =
reinterpret_cast< float*
> ( f );
delete foo; } );
222 pybind11::array_t < float > retArr = pybind11::array_t<float> ( {
static_cast<int> ( settings->
detectedSymmetry.size() ),
static_cast<int> ( 6 ) },
223 { 6 *
sizeof(float),
sizeof(
float) },
225 pyCapsuleStrRecSym );
229 },
"This function returns the recommended symmetry axes as a 2D numpy array." )
230 .def (
"getAllCSyms",
234 float* npVals =
new float[
static_cast<unsigned int> ( settings->
allDetectedCAxes.size() * 6 )];
238 for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( settings->
allDetectedCAxes.size() ); iter++ ) {
for ( proshade_unsign it = 0; it < 6; it++ ) { npVals[(iter*6)+it] = settings->
allDetectedCAxes.at(iter).at(it); } }
241 pybind11::capsule pyCapsuleStrSymList ( npVals, [](
void *f ) {
float* foo =
reinterpret_cast< float*
> ( f );
delete foo; } );
244 pybind11::array_t < float > retArr = pybind11::array_t<float> ( {
static_cast<int> ( settings->
allDetectedCAxes.size() ), 6 },
245 { 6 *
sizeof(float),
sizeof(
float) },
247 pyCapsuleStrSymList );
251 },
"This function returns all symmetry axes as a 2D numpy array." )
252 .def (
"getNonCSymmetryAxesIndices",
256 pybind11::dict retDict;
257 pybind11::list dList, tList, oList, iList;
260 for ( proshade_unsign dIt = 0; dIt < static_cast<proshade_unsign> ( settings->
allDetectedDAxes.size() ); dIt++ )
262 pybind11::list memList;
263 for ( proshade_unsign memIt = 0; memIt < static_cast<proshade_unsign> ( settings->
allDetectedDAxes.at(dIt).size() ); memIt++ )
267 dList.append ( memList );
271 for ( proshade_unsign tIt = 0; tIt < static_cast<proshade_unsign> ( settings->
allDetectedTAxes.size() ); tIt++ )
277 for ( proshade_unsign oIt = 0; oIt < static_cast<proshade_unsign> ( settings->
allDetectedOAxes.size() ); oIt++ )
283 for ( proshade_unsign iIt = 0; iIt < static_cast<proshade_unsign> ( settings->
allDetectedIAxes.size() ); iIt++ )
289 retDict[ pybind11::handle ( pybind11::str (
"D" ).ptr ( ) ) ] = dList;
290 retDict[ pybind11::handle ( pybind11::str (
"T" ).ptr ( ) ) ] = tList;
291 retDict[ pybind11::handle ( pybind11::str (
"O" ).ptr ( ) ) ] = oList;
292 retDict[ pybind11::handle ( pybind11::str (
"I" ).ptr ( ) ) ] = iList;
296 },
"This function returns array of non-C axes indices." )
297 .def (
"getAllGroupElements",
301 pybind11::buffer_info buf = axList.request();
302 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 ); }
305 std::vector< proshade_unsign > axesList;
309 std::vector < std::vector< proshade_double > > vals =
self.getAllGroupElements ( settings, axesList, groupType, matrixTolerance );
312 pybind11::list retList;
315 for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( vals.size() ); iter++ )
318 float* npVals =
new float[
static_cast<unsigned int> ( 9 )];
322 for ( proshade_unsign it = 0; it < 9; it++ ) { npVals[it] = vals.at(iter).at(it); }
325 pybind11::capsule pyCapsuleGrpEl ( npVals, [](
void *f ) {
float* foo =
reinterpret_cast< float*
> ( f );
delete foo; } );
328 pybind11::array_t < float > retArr = pybind11::array_t<float> ( { 3, 3 },
329 { 3 *
sizeof(float),
sizeof(
float) },
334 retList.append ( retArr );
339 },
"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 (
"getBestRotationMapPeaksEulerAngles",
347 std::vector< proshade_double > vals =
self.getBestRotationMapPeaksEulerAngles ( settings );
350 float* npVals =
new float[
static_cast<unsigned int> (vals.size())];
354 for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( vals.size() ); iter++ ) { npVals[iter] = vals.at(iter); }
357 pybind11::capsule pyCapsuleEuPeak ( npVals, [](
void *f ) {
float* foo =
reinterpret_cast< float*
> ( f );
delete foo; } );
360 pybind11::array_t < float > retArr = pybind11::array_t<float> ( {
static_cast<unsigned int> (vals.size()) },
367 },
"This function returns a vector of three floats, the three Euler angles of the best peak in the rotation map.", pybind11::arg (
"settings" ) )
368 .def (
"getBestRotationMapPeaksRotationMatrix",
372 std::vector< proshade_double > vals =
self.getBestRotationMapPeaksEulerAngles ( settings );
375 proshade_double* retMat =
new proshade_double[9];
380 pybind11::capsule pyCapsuleRMPeak ( retMat, [](
void *f ) { proshade_double* foo =
reinterpret_cast< proshade_double*
> ( f );
delete foo; } );
383 pybind11::array_t < proshade_double > retArr = pybind11::array_t<proshade_double> ( { 3, 3 },
384 { 3 *
sizeof(proshade_double),
sizeof(proshade_double) },
390 },
"This function returns a rotation matrix representing the best peak in the rotation map.", pybind11::arg (
"settings" ) )
391 .def (
"rotateMap", &
ProSHADE_internal_data::ProSHADE_data::rotateMap,
"This function rotates a map based on the given Euler angles.", pybind11::arg (
"settings" ), pybind11::arg (
"eulerAlpha" ), pybind11::arg (
"eulerBeta" ), pybind11::arg (
"eulerGamma" ) )
394 .def (
"getOverlayTranslations",
398 std::vector< proshade_double > vals =
self.getBestTranslationMapPeaksAngstrom ( staticStructure, eulA, eulB, eulG );
401 pybind11::dict retDict;
402 pybind11::list rotCen, toOverlay;
405 rotCen.append (
self.originalPdbRotCenX );
406 rotCen.append (
self.originalPdbRotCenY );
407 rotCen.append (
self.originalPdbRotCenZ );
409 toOverlay.append (
self.originalPdbTransX );
410 toOverlay.append (
self.originalPdbTransY );
411 toOverlay.append (
self.originalPdbTransZ );
414 retDict[ pybind11::handle ( pybind11::str (
"centreOfRotation" ).ptr ( ) ) ] = rotCen;
415 retDict[ pybind11::handle ( pybind11::str (
"rotCenToOverlay" ).ptr ( ) ) ] = toOverlay;
419 },
"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" ) )
420 .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 (
"settings" ), pybind11::arg (
"trsX" ), pybind11::arg (
"trsY" ), pybind11::arg (
"trsZ" ) )
423 .def (
"findSHIndex",
427 proshade_signed index = seanindex ( order, band,
self.spheres[shell]->getLocalBandwidth() );
431 },
"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." )
432 .def (
"getSphericalHarmonics",
436 std::complex<proshade_double>* npVals =
new std::complex<proshade_double>[
static_cast<proshade_unsign
> (
self.noSpheres * pow(
self.maxShellBand, 2.0 ) )];
440 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 ); }
443 proshade_signed pyPosSH;
444 proshade_signed pyPos;
447 for ( proshade_signed shIt = 0; shIt < static_cast<proshade_signed> (
self.noSpheres ); shIt++ )
449 for ( proshade_signed bnd = 0; bnd < static_cast<proshade_signed> (
self.spheres[shIt]->getLocalBandwidth() ); bnd++ )
451 for ( proshade_signed order = -bnd; order <= bnd; order++ )
453 pyPosSH = ( shIt * pow(
self.maxShellBand, 2.0 ) );
454 pyPos = seanindex ( order, bnd,
self.spheres[shIt]->getLocalBandwidth() );
455 npVals[pyPosSH+pyPos].real (
self.sphericalHarmonics[shIt][pyPos][0] );
456 npVals[pyPosSH+pyPos].imag (
self.sphericalHarmonics[shIt][pyPos][1] );
462 pybind11::capsule pyCapsuleSHs ( npVals, [](
void *f ) { std::complex<proshade_double>* foo =
reinterpret_cast< std::complex<proshade_double>*
> ( f );
delete foo; } );
465 pybind11::array_t < std::complex<proshade_double> > retArr = pybind11::array_t < std::complex<proshade_double > > (
466 {
static_cast<int> (
self.noSpheres ),
static_cast<int> ( pow (
self.maxShellBand, 2.0 ) ) },
467 {
sizeof ( std::complex < proshade_double > ) *
static_cast<int> ( pow (
self.maxShellBand, 2.0 ) ),
sizeof ( std::complex < proshade_double > ) },
473 },
"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." )
478 std::complex<proshade_double>* npVals =
new std::complex < proshade_double >[
static_cast<proshade_unsign
> (
self.maxShellBand * ( (
self.maxShellBand * 2 ) + 1 ) * ( (
self.maxShellBand * 2 ) + 1 ) )];
482 proshade_signed index = 0;
485 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 ); }
488 for ( proshade_signed bandIter = 0; bandIter < static_cast< proshade_signed > (
self.maxShellBand ); bandIter++ )
490 for ( proshade_unsign order1 = 0; order1 < ( ( bandIter * 2 ) + 1 ); order1++ )
492 for ( proshade_unsign order2 = 0; order2 < ( ( bandIter * 2 ) + 1 ); order2++ )
494 index = order2 + ( (
self.maxShellBand * 2 ) + 1 ) * ( order1 + ( (
self.maxShellBand * 2 ) + 1 ) * bandIter );
495 npVals[index].real (
self.eMatrices[bandIter][order1][order2][0] );
496 npVals[index].imag (
self.eMatrices[bandIter][order1][order2][1] );
502 pybind11::capsule pyCapsuleEMs ( npVals, [](
void *f ) { std::complex<proshade_double>* foo =
reinterpret_cast< std::complex<proshade_double>*
> ( f );
delete foo; } );
505 pybind11::array_t < std::complex<proshade_double> > retArr = pybind11::array_t < std::complex<proshade_double > > (
506 {
self.maxShellBand, ( (
self.maxShellBand * 2 ) + 1 ), ( (
self.maxShellBand * 2 ) + 1 ) },
507 {
sizeof ( std::complex < proshade_double > ) *
static_cast<int> ( ( (
self.maxShellBand * 2 ) + 1 ) * ( (
self.maxShellBand * 2 ) + 1 ) ),
508 sizeof ( std::complex < proshade_double > ) *
static_cast<int> ( (
self.maxShellBand * 2 ) + 1 ),
509 sizeof ( std::complex < proshade_double > ) },
515 },
"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." )
516 .def (
"getSO3Coefficients",
520 std::complex<proshade_double>* npVals =
new std::complex < proshade_double >[
static_cast<proshade_unsign
> (
self.maxShellBand * ( (
self.maxShellBand * 2 ) + 1 ) * ( (
self.maxShellBand * 2 ) + 1 ) )];
524 proshade_signed index = 0;
527 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 ); }
530 for ( proshade_signed bandIter = 0; bandIter < static_cast< proshade_signed > (
self.maxShellBand ); bandIter++ )
532 for ( proshade_unsign order1 = 0; order1 < ( ( bandIter * 2 ) + 1 ); order1++ )
534 for ( proshade_unsign order2 = 0; order2 < ( ( bandIter * 2 ) + 1 ); order2++ )
536 index = order2 + ( (
self.maxShellBand * 2 ) + 1 ) * ( order1 + ( (
self.maxShellBand * 2 ) + 1 ) * bandIter );
537 npVals[index].real (
self.so3Coeffs[
self.so3CoeffsArrayIndex ( order1 - bandIter, order2 - bandIter, bandIter )][0] );
538 npVals[index].imag (
self.so3Coeffs[
self.so3CoeffsArrayIndex ( order1 - bandIter, order2 - bandIter, bandIter )][1] );
544 pybind11::capsule pyCapsuleSOCoeffs ( npVals, [](
void *f ) { std::complex<proshade_double>* foo =
reinterpret_cast< std::complex<proshade_double>*
> ( f );
delete foo; } );
547 pybind11::array_t < std::complex<proshade_double> > retArr = pybind11::array_t < std::complex<proshade_double > > (
548 {
self.maxShellBand, ( (
self.maxShellBand * 2 ) + 1 ), ( (
self.maxShellBand * 2 ) + 1 ) },
549 {
sizeof ( std::complex < proshade_double > ) *
static_cast<int> ( ( (
self.maxShellBand * 2 ) + 1 ) * ( (
self.maxShellBand * 2 ) + 1 ) ),
550 sizeof ( std::complex < proshade_double > ) *
static_cast<int> ( (
self.maxShellBand * 2 ) + 1 ),
551 sizeof ( std::complex < proshade_double > ) },
557 },
"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." )
558 .def (
"getRotationFunctionMap",
562 std::complex<proshade_double>* npVals =
new std::complex < proshade_double >[
static_cast<proshade_unsign
> ( (
self.maxShellBand * 2 ) * (
self.maxShellBand * 2 ) * (
self.maxShellBand * 2 ) )];
566 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] ); }
569 pybind11::capsule pyCapsuleRotMap ( npVals, [](
void *f ) { std::complex<proshade_double>* foo =
reinterpret_cast< std::complex<proshade_double>*
> ( f );
delete foo; } );
572 pybind11::array_t < std::complex < proshade_double > > retArr = pybind11::array_t < std::complex < proshade_double > > (
573 { (
self.maxShellBand * 2 ), (
self.maxShellBand * 2 ), (
self.maxShellBand * 2 ) },
574 { (
self.maxShellBand * 2 ) * (
self.maxShellBand * 2 ) *
sizeof(std::complex < proshade_double >),
575 (
self.maxShellBand * 2 ) *
sizeof(std::complex < proshade_double >),
576 sizeof(std::complex < proshade_double >) },
582 },
"This function returns the (self) rotation function as a three-dimensional map of complex numbers." )
583 .def (
"getRotationMatrixFromSOFTCoordinates",
587 proshade_double* npVals =
new proshade_double[9];
591 proshade_double eulA, eulB, eulG;
600 pybind11::capsule pyCapsuleRMSoft ( npVals, [](
void *f ) { proshade_double* foo =
reinterpret_cast< proshade_double*
> ( f );
delete foo; } );
603 pybind11::array_t < proshade_double > retArr = pybind11::array_t < proshade_double > ( { 3, 3 },
604 { 3 *
sizeof(proshade_double),
sizeof(proshade_double) },
610 },
"This function converts a given rotation function map position onto the corresponding rotation matrix.", pybind11::arg (
"xPos" ), pybind11::arg (
"yPos" ), pybind11::arg (
"zPos" ) )
611 .def (
"getTranslationFunctionMap",
615 std::complex<proshade_double>* npVals =
new std::complex < proshade_double >[
static_cast<proshade_unsign
> (
self.getXDim() *
self.getYDim() *
self.getZDim() )];
619 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] ); }
622 pybind11::capsule pyCapsuleTrsMap ( npVals, [](
void *f ) { std::complex<proshade_double>* foo =
reinterpret_cast< std::complex<proshade_double>*
> ( f );
delete foo; } );
625 pybind11::array_t < std::complex < proshade_double > > retArr = pybind11::array_t < std::complex < proshade_double > > (
626 {
self.getXDim(),
self.getYDim(),
self.getZDim() },
627 {
self.getYDim() *
self.getZDim() *
sizeof(std::complex < proshade_double >),
628 self.getZDim() *
sizeof(std::complex < proshade_double >),
629 sizeof(std::complex < proshade_double >) },
635 },
"This function returns the translation function as a three-dimensional map of complex numbers." )
673 .def (
"__repr__", [] (
const ProSHADE_internal_data::ProSHADE_data &a ) {
return "<ProSHADE_data class object> (This class contains all information, results and available functionalities for a structure)"; } );
676 pyProSHADE.def (
"computeGroupElementsForGroup",
677 [] ( proshade_double x, proshade_double y, proshade_double z, proshade_unsign fold ) -> pybind11::array_t < proshade_double >
683 proshade_double* npVals =
new proshade_double[
static_cast<proshade_unsign
> ( retVec.size() ) * 9];
687 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); } }
690 pybind11::capsule pyCapsuleGrEls ( npVals, [](
void *f ) { proshade_double* foo =
reinterpret_cast< proshade_double*
> ( f );
delete foo; } );
693 pybind11::array_t < proshade_double > retArr = pybind11::array_t < proshade_double > (
694 {
static_cast< int > ( retVec.size() ), 3, 3 },
695 { 9 *
sizeof(proshade_double), 3 *
sizeof(proshade_double),
sizeof(proshade_double) },
701 },
"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" ) );
702 pyProSHADE.def (
"joinElementsFromDifferentGroups",
703 [] ( pybind11::array_t < proshade_double, pybind11::array::c_style | pybind11::array::forcecast > first,
704 pybind11::array_t < proshade_double, pybind11::array::c_style | pybind11::array::forcecast > second,
705 proshade_double matrixTolerance,
706 bool combine ) -> pybind11::array_t < proshade_double >
709 pybind11::buffer_info buf1 = first.request();
710 pybind11::buffer_info buf2 = second.request();
713 if ( buf1.ndim != 3 || buf2.ndim != 3 )
715 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;
716 return ( pybind11::array_t < proshade_double > () );
720 std::vector< std::vector< proshade_double > > fVec, sVec;
722 proshade_double* dataPtr1 =
reinterpret_cast < proshade_double*
> ( buf1.ptr );
723 for ( proshade_unsign elIt = 0; elIt < static_cast<proshade_unsign> ( buf1.shape.at(0) ); elIt++ )
725 std::vector< proshade_double > rotMat;
726 for ( proshade_unsign rowIt = 0; rowIt < static_cast<proshade_unsign> ( buf1.shape.at(1) ); rowIt++ )
728 for ( proshade_unsign colIt = 0; colIt < static_cast<proshade_unsign> ( buf1.shape.at(2) ); colIt++ )
736 proshade_double* dataPtr2 =
reinterpret_cast < proshade_double*
> ( buf2.ptr );
737 for ( proshade_unsign elIt = 0; elIt < static_cast<proshade_unsign> ( buf2.shape.at(0) ); elIt++ )
739 std::vector< proshade_double > rotMat;
740 for ( proshade_unsign rowIt = 0; rowIt < static_cast<proshade_unsign> ( buf2.shape.at(1) ); rowIt++ )
742 for ( proshade_unsign colIt = 0; colIt < static_cast<proshade_unsign> ( buf2.shape.at(2) ); colIt++ )
755 proshade_double* npVals =
new proshade_double[
static_cast<proshade_unsign
> ( retVec.size() ) * 9];
759 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); } }
762 pybind11::capsule pyCapsuleGrElsCombo ( npVals, [](
void *f ) { proshade_double* foo =
reinterpret_cast< proshade_double*
> ( f );
delete foo; } );
765 pybind11::array_t < proshade_double > retArr = pybind11::array_t < proshade_double > (
766 {
static_cast< int > ( retVec.size() ), 3, 3 },
767 { 9 *
sizeof(proshade_double), 3 *
sizeof(proshade_double),
sizeof(proshade_double) },
769 pyCapsuleGrElsCombo );
773 },
"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" ) );