ProSHADE  0.7.6.2 (DEC 2021)
Protein Shape Detection
ProSHADE_internal_io Namespace Reference

This namespace contains the internal input/output functions. None of these should be directly accessed by the user. More...

Enumerations

enum  InputType { UNKNOWN, PDB, MAP, GEMMI }
 

Functions

InputType figureDataType (std::string fName)
 Function determining input data type. More...
 
bool isFilePDB (std::string fName)
 Function determining if the input data type is PDB. More...
 
bool isFileMAP (std::string fName)
 Function determining if the input data type is MAP. More...
 
void readInMapHeader (gemmi::Ccp4< float > *map, proshade_unsign *xDimInds, proshade_unsign *yDimInds, proshade_unsign *zDimInds, proshade_single *xDim, proshade_single *yDim, proshade_single *zDim, proshade_single *aAng, proshade_single *bAng, proshade_single *cAng, proshade_signed *xFrom, proshade_signed *yFrom, proshade_signed *zFrom, proshade_signed *xAxOrigin, proshade_signed *yAxOrigin, proshade_signed *zAxOrigin, proshade_unsign *xAxOrder, proshade_unsign *yAxOrder, proshade_unsign *zAxOrder, proshade_unsign *xGridInds, proshade_unsign *yGridInds, proshade_unsign *zGridInds)
 This function parses the CCP4 MAP file header as read in by gemmi. More...
 
void readInMapData (gemmi::Ccp4< float > *gemmiMap, proshade_double *&map, proshade_unsign xDimInds, proshade_unsign yDimInds, proshade_unsign zDimInds, proshade_unsign xAxOrder, proshade_unsign yAxOrder, proshade_unsign zAxOrder)
 This function converts the gemmi Ccp4 object data to ProSHADE internal map representation. More...
 
void readInMapData (gemmi::Ccp4< int8_t > *gemmiMap, proshade_double *&map, proshade_unsign xDimInds, proshade_unsign yDimInds, proshade_unsign zDimInds, proshade_unsign xAxOrder, proshade_unsign yAxOrder, proshade_unsign zAxOrder)
 This function converts the gemmi Ccp4 object data to ProSHADE mask representation. More...
 
void applyMask (proshade_double *&map, std::string maskFile, proshade_unsign xDimInds, proshade_unsign yDimInds, proshade_unsign zDimInds, proshade_signed verbose, proshade_signed messageShift, proshade_double *maskArray=nullptr, proshade_unsign maXInds=0, proshade_unsign maYInds=0, proshade_unsign maZInds=0)
 This function reads and applies the mask to the map. More...
 
void applyMaskFromArray (proshade_double *&map, proshade_unsign xDimInds, proshade_unsign yDimInds, proshade_unsign zDimInds, proshade_double *&mask, proshade_unsign xDimIndsMsk, proshade_unsign yDimIndsMsk, proshade_unsign zDimIndsMsk, proshade_signed verbose, proshade_signed messageShift)
 This function applies the mask to the map. More...
 
void applyWeights (proshade_double *&map, std::string weightsFile, proshade_unsign xDimInds, proshade_unsign yDimInds, proshade_unsign zDimInds, proshade_signed verbose, proshade_signed messageShift, proshade_double *weightsArray=nullptr, proshade_unsign waXInds=0, proshade_unsign waYInds=0, proshade_unsign waZInds=0)
 This function reads and applies the Fourier weights to the map. More...
 
void applyWeightsFromArray (proshade_double *&map, proshade_unsign xDimInds, proshade_unsign yDimInds, proshade_unsign zDimInds, proshade_double *&weights, proshade_unsign xDimIndsWgh, proshade_unsign yDimIndsWgh, proshade_unsign zDimIndsWgh, proshade_signed verbose, proshade_signed messageShift)
 This function applies the weights to the map. More...
 
void writeOutMapHeader (gemmi::Ccp4< float > *map, proshade_unsign xDimInds, proshade_unsign yDimInds, proshade_unsign zDimInds, proshade_single xDim, proshade_single yDim, proshade_single zDim, proshade_single aAng, proshade_single bAng, proshade_single cAng, proshade_signed xFrom, proshade_signed yFrom, proshade_signed zFrom, proshade_signed xAxOrigin, proshade_signed yAxOrigin, proshade_signed zAxOrigin, proshade_unsign xAxOrder, proshade_unsign yAxOrder, proshade_unsign zAxOrder, proshade_unsign xGridInds, proshade_unsign yGridInds, proshade_unsign zGridInds, std::string title, int mode)
 This function parses the CCP4 MAP file header as read in by gemmi. More...
 
void writeRotationTranslationJSON (proshade_double trsX1, proshade_double trsY1, proshade_double trsZ1, proshade_double eulA, proshade_double eulB, proshade_double eulG, proshade_double trsX2, proshade_double trsY2, proshade_double trsZ2, std::string fileName)
 Function for writing out the optimal rotation and translation into a JSON file. More...
 

Detailed Description

This namespace contains the internal input/output functions. None of these should be directly accessed by the user.

The ProSHADE_internal_io namespace contains the helper functions for the data input and output. These should never be directly used by the user and these only serve to allow for self-documenting nature of the code. They are called internally by more advanced functions from the higher complexity classes.

Function Documentation

◆ applyMask()

void ProSHADE_internal_io::applyMask ( proshade_double *&  map,
std::string  maskFile,
proshade_unsign  xDimInds,
proshade_unsign  yDimInds,
proshade_unsign  zDimInds,
proshade_signed  verbose,
proshade_signed  messageShift,
proshade_double *  maskArray = nullptr,
proshade_unsign  maXInds = 0,
proshade_unsign  maYInds = 0,
proshade_unsign  maZInds = 0 
)

This function reads and applies the mask to the map.

This function reads in the mask file and determines, if its dimensions are the same as the input map or not. If not, then the mask is re-sampled is such a way as to have the same number of indices as the input map. Once the map and mask dimensions are the same, the mask is applied by multiplying all map values by corresponding mask value. The map change is done in-place!

Parameters
[in]mapPointer to a the internal variable containing the map (it will be changed in place).
[in]maskFileFilename of where the mask is to be read from.
[in]xDimIndsThe size of x dimension in indices.
[in]yDimIndsThe size of y dimension in indices.
[in]zDimIndsThe size of z dimension in indices.
[in]verboseHow much std::out output would you like?
[in]messageShiftAre we in a subprocess, so that the log should be shifted for this function call? If so, by how much?
[in]maskArrayAn array of mask values (default nullptr) to be used instead of an input file.
[in]maXIndsThe size of maskArray x dimension in indices (defaults to 0).
[in]maYIndsThe size of maskArray y dimension in indices (defaults to 0).
[in]maZIndsThe size of maskArray z dimension in indices (defaults to 0).

Definition at line 295 of file ProSHADE_io.cpp.

296 {
297  //================================================ Report progress
298  std::stringstream hlpSS;
299  hlpSS << "Reading mask " << maskFile;
300  ProSHADE_internal_messages::printProgressMessage ( verbose, 2, hlpSS.str() );
301 
302  //================================================ Are we reading from array or from file?
303  if ( ( maskArray != nullptr ) && ( maXInds != 0 ) && ( maYInds != 0 ) && ( maZInds != 0 ) )
304  {
305  //============================================ Array it is!
306  ProSHADE_internal_io::applyMaskFromArray ( map, xDimInds, yDimInds, zDimInds, maskArray, maXInds, maYInds, maZInds, verbose, messageShift );
307  }
308  else
309  {
310  //============================================ Check if filename was given
311  if ( maskFile == "" ) { ProSHADE_internal_messages::printProgressMessage ( verbose, 3, "No mask supplied." ); return; }
312 
313  //============================================ File it is! Open the mask
314  gemmi::Ccp4<float> mask;
315  mask.read_ccp4 ( gemmi::MaybeGzipped ( maskFile.c_str() ) );
316 
317  //============================================ Convert to XYZ and create complete mask, if need be
318  mask.setup ( gemmi::GridSetup::ReorderOnly, 0 );
319 
320  //============================================ Read in the rest of the mask file header
321  proshade_unsign xDI, yDI, zDI, xAOR, yAOR, zAOR, xGI, yGI, zGI;
322  proshade_single xDS, yDS, zDS, aA, bA, cA;
323  proshade_signed xF, yF, zF, xAO, yAO, zAO;
325  &xDI, &yDI, &zDI,
326  &xDS, &yDS, &zDS,
327  &aA, &bA, &cA,
328  &xF, &yF, &zF,
329  &xAO, &yAO, &zAO,
330  &xAOR, &yAOR, &zAOR,
331  &xGI, &yGI, &zGI );
332 
333  //============================================ Save the mask values to ProSHADE variable
334  proshade_double* internalMask = nullptr;
335  ProSHADE_internal_io::readInMapData ( &mask, internalMask, xDI, yDI, zDI, xAOR, yAOR, zAOR );
336 
337  //============================================ Apply mask from array
338  ProSHADE_internal_io::applyMaskFromArray ( map, xDimInds, yDimInds, zDimInds, internalMask, xDI, yDI, zDI, verbose, messageShift );
339 
340  //============================================ Release the memory
341  delete[] internalMask;
342  }
343 
344  //================================================ Done
345  return ;
346 
347 }

◆ applyMaskFromArray()

void ProSHADE_internal_io::applyMaskFromArray ( proshade_double *&  map,
proshade_unsign  xDimInds,
proshade_unsign  yDimInds,
proshade_unsign  zDimInds,
proshade_double *&  mask,
proshade_unsign  xDimIndsMsk,
proshade_unsign  yDimIndsMsk,
proshade_unsign  zDimIndsMsk,
proshade_signed  verbose,
proshade_signed  messageShift 
)

This function applies the mask to the map.

This function firstlydetermines if its dimensions of mask are the same as the input map or not. If not, then the mask is re-sampled is such a way as to have the same number of indices as the input map. Once the map and mask dimensions are the same, the mask is applied by multiplying all map values by corresponding mask value. The map change is done in-place!

Parameters
[in]mapPointer to a the internal variable containing the map (it will be changed in place).
[in]xDimIndsThe size of the map x dimension in indices.
[in]yDimIndsThe size of the map y dimension in indices.
[in]zDimIndsThe size of the map z dimension in indices.
[in]maskA pointer to 1D array of the mask values.
[in]xDimIndsMskThe size of the mask x dimension in indices.
[in]yDimIndsMskThe size of the mask y dimension in indices.
[in]zDimIndsMskThe size of the mask z dimension in indices.
[in]verboseHow much std::out output would you like?
[in]messageShiftAre we in a subprocess, so that the log should be shifted for this function call? If so, by how much?

Definition at line 367 of file ProSHADE_io.cpp.

368 {
369  //================================================ Initialise local variables
370  size_t origVolume = xDimInds * yDimInds * zDimInds;
371  size_t newVolume = xDimIndsMsk * yDimIndsMsk * zDimIndsMsk;
372  proshade_double* maskFinal;
373 
374  //================================================ If mask has different number of indices than map, then re-sample mask
375  if ( ( xDimIndsMsk != xDimInds ) || ( yDimIndsMsk != yDimInds ) || ( zDimIndsMsk != zDimInds ) )
376  {
377  //============================================ Initialise variables
378  fftw_complex* origCoeffs = new fftw_complex [newVolume ];
379  fftw_complex* origCoeffsHKL = new fftw_complex [newVolume ];
380  fftw_complex* modifCoeffs = new fftw_complex [origVolume];
381  fftw_complex* modifCoeffsHKL = new fftw_complex [origVolume];
382  fftw_complex* inMap = new fftw_complex [newVolume ];
383  fftw_complex* outMap = new fftw_complex [origVolume];
384 
385  //============================================ Check memory allocation
386  ProSHADE_internal_misc::checkMemoryAllocation ( inMap, __FILE__, __LINE__, __func__ );
387  ProSHADE_internal_misc::checkMemoryAllocation ( outMap, __FILE__, __LINE__, __func__ );
388  ProSHADE_internal_misc::checkMemoryAllocation ( origCoeffs, __FILE__, __LINE__, __func__ );
389  ProSHADE_internal_misc::checkMemoryAllocation ( modifCoeffs, __FILE__, __LINE__, __func__ );
390  ProSHADE_internal_misc::checkMemoryAllocation ( origCoeffsHKL, __FILE__, __LINE__, __func__ );
391  ProSHADE_internal_misc::checkMemoryAllocation ( modifCoeffsHKL, __FILE__, __LINE__, __func__ );
392 
393  //============================================ Set array to zeroes
394  for ( size_t iter = 0; iter < origVolume; iter++ ) { modifCoeffsHKL[iter][0] = 0.0; modifCoeffsHKL[iter][1] = 0.0; }
395 
396  //============================================ Cope mask to Fourier input array
397  for ( size_t iter = 0; iter < newVolume; iter++ ) { inMap[iter][0] = mask[iter]; inMap[iter][1] = 0.0; }
398 
399  //============================================ Prepare Fourier transform plans
400  fftw_plan planForwardFourier = fftw_plan_dft_3d ( static_cast< int > ( xDimIndsMsk ), static_cast< int > ( yDimIndsMsk ), static_cast< int > ( zDimIndsMsk ), inMap, origCoeffs, FFTW_FORWARD, FFTW_ESTIMATE );
401  fftw_plan inverseFoourier = fftw_plan_dft_3d ( static_cast< int > ( xDimInds ), static_cast< int > ( yDimInds ), static_cast< int > ( zDimInds ), modifCoeffs, outMap, FFTW_BACKWARD, FFTW_ESTIMATE );
402 
403  //============================================ Compute pre and post changes
404  proshade_signed xPre, yPre, zPre;
405  xPre = std::abs ( ( static_cast< proshade_signed > ( xDimInds ) - static_cast< proshade_signed > ( xDimIndsMsk ) ) / 2 );
406  yPre = std::abs ( ( static_cast< proshade_signed > ( yDimInds ) - static_cast< proshade_signed > ( yDimIndsMsk ) ) / 2 );
407  zPre = std::abs ( ( static_cast< proshade_signed > ( zDimInds ) - static_cast< proshade_signed > ( zDimIndsMsk ) ) / 2 );
408 
409  if ( ( ( static_cast< proshade_signed > ( xDimInds ) - static_cast< proshade_signed > ( xDimIndsMsk ) ) % 2 ) == 1 ) { xPre -= 1; }
410  if ( ( ( static_cast< proshade_signed > ( yDimInds ) - static_cast< proshade_signed > ( yDimIndsMsk ) ) % 2 ) == 1 ) { yPre -= 1; }
411  if ( ( ( static_cast< proshade_signed > ( zDimInds ) - static_cast< proshade_signed > ( zDimIndsMsk ) ) % 2 ) == 1 ) { zPre -= 1; }
412 
413  //============================================ Run forward Fourier
414  fftw_execute ( planForwardFourier );
415 
416  //============================================ Initialise local variables
417  proshade_signed maskMapIndex = 0;
418  proshade_signed densMapIndex = 0;
419  proshade_signed xMaskPos, yMaskPos, zMaskPos, xDensPos, yDensPos, zDensPos;
420  proshade_signed maskH, maskK, maskL;
421 
422  //============================================ Convert mask to HKL for re-boxing
423  for ( proshade_signed xIt = 0; xIt < static_cast< proshade_signed > ( xDimIndsMsk ); xIt++ )
424  {
425  for ( proshade_signed yIt = 0; yIt < static_cast< proshade_signed > ( yDimIndsMsk ); yIt++ )
426  {
427  for ( proshade_signed zIt = 0; zIt < static_cast< proshade_signed > ( zDimIndsMsk ); zIt++ )
428  {
429  //================================ Convert to HKL
430  maskH = xIt + static_cast< proshade_signed > ( (xDimIndsMsk+1) / 2 ); if ( maskH >= static_cast< proshade_signed > ( xDimIndsMsk ) ) { maskH -= xDimIndsMsk; }
431  maskK = yIt + static_cast< proshade_signed > ( (yDimIndsMsk+1) / 2 ); if ( maskK >= static_cast< proshade_signed > ( yDimIndsMsk ) ) { maskK -= yDimIndsMsk; }
432  maskL = zIt + static_cast< proshade_signed > ( (zDimIndsMsk+1) / 2 ); if ( maskL >= static_cast< proshade_signed > ( zDimIndsMsk ) ) { maskL -= zDimIndsMsk; }
433 
434  //================================ Find the positions
435  maskMapIndex = zIt + static_cast< proshade_signed > ( zDimIndsMsk ) * ( yIt + static_cast< proshade_signed > ( yDimIndsMsk ) * xIt );
436  densMapIndex = maskL + static_cast< proshade_signed > ( zDimIndsMsk ) * ( maskK + static_cast< proshade_signed > ( yDimIndsMsk ) * maskH );
437 
438  //================================ Save the values
439  origCoeffsHKL[densMapIndex][0] = origCoeffs[maskMapIndex][0];
440  origCoeffsHKL[densMapIndex][1] = origCoeffs[maskMapIndex][1];
441  }
442  }
443  }
444 
445  //============================================ Rebox
446  for ( proshade_signed xIt = 0; xIt < static_cast< proshade_signed > ( xDimInds ); xIt++ )
447  {
448  for ( proshade_signed yIt = 0; yIt < static_cast< proshade_signed > ( yDimInds ); yIt++ )
449  {
450  for ( proshade_signed zIt = 0; zIt < static_cast< proshade_signed > ( zDimInds ); zIt++ )
451  {
452  //================================ Deal with X
453  if ( xDimIndsMsk >= xDimInds ) { xMaskPos = xIt + xPre; }
454  else { xMaskPos = xIt - xPre; }
455  xDensPos = xIt;
456 
457  //================================ Deal with Y
458  if ( yDimIndsMsk >= yDimInds ) { yMaskPos = yIt + yPre; }
459  else { yMaskPos = yIt - yPre; }
460  yDensPos = yIt;
461 
462  //================================ Deal with Z
463  if ( zDimIndsMsk >= zDimInds ) { zMaskPos = zIt + zPre; }
464  else { zMaskPos = zIt - zPre; }
465  zDensPos = zIt;
466 
467  //================================ Skip if mask value not available (because the modifCoeffsHKL array is zeroed, we do not need to do anything here)
468  if ( ( xMaskPos < 0 ) || ( xMaskPos >= static_cast< proshade_signed > ( xDimIndsMsk ) ) ) { continue; }
469  if ( ( yMaskPos < 0 ) || ( yMaskPos >= static_cast< proshade_signed > ( yDimIndsMsk ) ) ) { continue; }
470  if ( ( zMaskPos < 0 ) || ( zMaskPos >= static_cast< proshade_signed > ( zDimIndsMsk ) ) ) { continue; }
471 
472  //================================ Find the positions
473  maskMapIndex = zMaskPos + static_cast< proshade_signed > ( zDimIndsMsk ) * ( yMaskPos + static_cast< proshade_signed > ( yDimIndsMsk ) * xMaskPos );
474  densMapIndex = zDensPos + static_cast< proshade_signed > ( zDimInds ) * ( yDensPos + static_cast< proshade_signed > ( yDimInds ) * xDensPos );
475 
476  //================================ Copy values
477  modifCoeffsHKL[densMapIndex][0] = origCoeffsHKL[maskMapIndex][0];
478  modifCoeffsHKL[densMapIndex][1] = origCoeffsHKL[maskMapIndex][1];
479  }
480  }
481  }
482 
483  //============================================ Convert mask back to FFTW order
484  for ( proshade_signed xIt = 0; xIt < static_cast< proshade_signed > ( xDimInds ); xIt++ )
485  {
486  for ( proshade_signed yIt = 0; yIt < static_cast< proshade_signed > ( yDimInds ); yIt++ )
487  {
488  for ( proshade_signed zIt = 0; zIt < static_cast< proshade_signed > ( zDimInds ); zIt++ )
489  {
490  //================================ Convert to HKL
491  maskH = xIt + static_cast< proshade_signed > ( xDimInds / 2 ); if ( maskH >= static_cast< proshade_signed > ( xDimInds ) ) { maskH -= xDimInds; }
492  maskK = yIt + static_cast< proshade_signed > ( yDimInds / 2 ); if ( maskK >= static_cast< proshade_signed > ( yDimInds ) ) { maskK -= yDimInds; }
493  maskL = zIt + static_cast< proshade_signed > ( zDimInds / 2 ); if ( maskL >= static_cast< proshade_signed > ( zDimInds ) ) { maskL -= zDimInds; }
494 
495  //================================ Find the positions
496  maskMapIndex = zIt + static_cast< proshade_signed > ( zDimInds ) * ( yIt + static_cast< proshade_signed > ( yDimInds ) * xIt );
497  densMapIndex = maskL + static_cast< proshade_signed > ( zDimInds ) * ( maskK + static_cast< proshade_signed > ( yDimInds ) * maskH );
498 
499  //================================ Save the values
500  modifCoeffs[densMapIndex][0] = modifCoeffsHKL[maskMapIndex][0];
501  modifCoeffs[densMapIndex][1] = modifCoeffsHKL[maskMapIndex][1];
502  }
503  }
504  }
505 
506  //============================================ Run inverse Fourier on the modified coefficients
507  fftw_execute ( inverseFoourier );
508 
509  //============================================ Delete old mask and allocate memory for the new, re-sampled mask
510  maskFinal = new proshade_double [origVolume];
511  ProSHADE_internal_misc::checkMemoryAllocation ( maskFinal, __FILE__, __LINE__, __func__ );
512 
513  //======================================== Copy results into a new, properly sampled mask
514  for ( size_t iter = 0; iter < origVolume; iter++ ) { maskFinal[iter] = outMap[iter][0]; }
515 
516  //======================================== Release remaining memory
517  fftw_destroy_plan ( planForwardFourier );
518  fftw_destroy_plan ( inverseFoourier );
519  delete[] origCoeffs;
520  delete[] modifCoeffs;
521  delete[] origCoeffsHKL;
522  delete[] modifCoeffsHKL;
523  delete[] inMap;
524  delete[] outMap;
525  }
526  else
527  {
528  maskFinal = new proshade_double [origVolume];
529  ProSHADE_internal_misc::checkMemoryAllocation ( maskFinal, __FILE__, __LINE__, __func__ );
530  for ( size_t iter = 0; iter < origVolume; iter++ ) { maskFinal[iter] = mask[iter]; }
531  }
532 
533  //================================================ Apply the mask to the map
534  for ( size_t iter = 0; iter < static_cast< size_t > ( xDimInds * yDimInds * zDimInds ); iter++ ) { map[iter] *= maskFinal[iter]; }
535 
536  //================================================ Release memory
537  delete[] maskFinal;
538 
539  //================================================ Report progress
540  ProSHADE_internal_messages::printProgressMessage ( verbose, 3, "Mask read in and applied successfully.", messageShift );
541 
542  //================================================ Done
543  return ;
544 
545 }

◆ applyWeights()

void ProSHADE_internal_io::applyWeights ( proshade_double *&  map,
std::string  weightsFile,
proshade_unsign  xDimInds,
proshade_unsign  yDimInds,
proshade_unsign  zDimInds,
proshade_signed  verbose,
proshade_signed  messageShift,
proshade_double *  weightsArray = nullptr,
proshade_unsign  waXInds = 0,
proshade_unsign  waYInds = 0,
proshade_unsign  waZInds = 0 
)

This function reads and applies the Fourier weights to the map.

This function reads in the weights file and determines, if its dimensions are the same as the input map or not. If not, then the weights are re-sampled is such a way as to have the same number of indices as the input map. Once the map and weights map dimensions are the same, the weights are applied by multiplying all map values by corresponding weight value. The map change is done in-place!

Parameters
[in]mapPointer to a the internal variable containing the map (it will be changed in place).
[in]weightsFileFilename of where the mask is to be read from.
[in]xDimIndsThe size of x dimension in indices.
[in]yDimIndsThe size of y dimension in indices.
[in]zDimIndsThe size of z dimension in indices.
[in]verboseHow much std::out output would you like?
[in]messageShiftAre we in a subprocess, so that the log should be shifted for this function call? If so, by how much?
[in]weightsArrayAn array of weights (default nullptr) to be used instead of input file.
[in]waXIndsThe size of weightsArray x dimension in indices (defaults to 0).
[in]waYIndsThe size of weightsArray y dimension in indices (defaults to 0).
[in]waZIndsThe size of weightsArray z dimension in indices (defaults to 0).

Definition at line 566 of file ProSHADE_io.cpp.

567 {
568  //================================================ Report progress
569  std::stringstream hlpSS;
570  hlpSS << "Reading weights " << weightsFile;
571  ProSHADE_internal_messages::printProgressMessage ( verbose, 2, hlpSS.str() );
572 
573  //================================================ Are we reading from file, or from array?
574  if ( ( weightsArray != nullptr ) && ( waXInds != 0 ) && ( waYInds != 0 ) && ( waZInds != 0 ) )
575  {
576  //============================================ From array it is!
577  ProSHADE_internal_io::applyWeightsFromArray ( map, xDimInds, yDimInds, zDimInds, weightsArray, waXInds, waYInds, waZInds, verbose, messageShift );
578  }
579  else
580  {
581  //============================================ Check if weights file was given
582  if ( weightsFile == "" ) { ProSHADE_internal_messages::printProgressMessage ( verbose, 3, "No weights supplied. Assuming all weights to be 1.0.", messageShift ); return; }
583 
584  //============================================ From file it is! Open the weights file
585  gemmi::Ccp4<float> weights;
586  weights.read_ccp4 ( gemmi::MaybeGzipped ( weightsFile.c_str() ) );
587 
588  //============================================ Convert to XYZ and create complete weights, if need be
589  weights.setup ( gemmi::GridSetup::ReorderOnly, 0 );
590 
591  //============================================ Read in the rest of the weights file header
592  proshade_unsign xDI, yDI, zDI, xAOR, yAOR, zAOR, xGI, yGI, zGI;
593  proshade_single xDS, yDS, zDS, aA, bA, cA;
594  proshade_signed xF, yF, zF, xAO, yAO, zAO;
596  &xDI, &yDI, &zDI,
597  &xDS, &yDS, &zDS,
598  &aA, &bA, &cA,
599  &xF, &yF, &zF,
600  &xAO, &yAO, &zAO,
601  &xAOR, &yAOR, &zAOR,
602  &xGI, &yGI, &zGI );
603 
604  //============================================ Save the weights values to local variable
605  proshade_double* internalWeights = nullptr;
606  ProSHADE_internal_io::readInMapData ( &weights, internalWeights, xDI, yDI, zDI, xAOR, yAOR, zAOR );
607 
608  //============================================ Apply weights from array
609  ProSHADE_internal_io::applyWeightsFromArray ( map, xDimInds, yDimInds, zDimInds, internalWeights, xDI, yDI, zDI, verbose, messageShift );
610 
611  //============================================ Release the memory
612  delete[] internalWeights;
613  }
614 
615  //================================================ Done
616  return ;
617 
618 }

◆ applyWeightsFromArray()

void ProSHADE_internal_io::applyWeightsFromArray ( proshade_double *&  map,
proshade_unsign  xDimInds,
proshade_unsign  yDimInds,
proshade_unsign  zDimInds,
proshade_double *&  weights,
proshade_unsign  xDimIndsWgh,
proshade_unsign  yDimIndsWgh,
proshade_unsign  zDimIndsWgh,
proshade_signed  verbose,
proshade_signed  messageShift 
)

This function applies the weights to the map.

This function firstly determines if its dimensions of weights are the same as the input map or not. If not, then the weights are re-sampled is such a way as to have the same number of indices as the input map. Once the map and weights dimensions are the same, the weights are applied by multiplying all map Fourier coefficients by corresponding weights value. The map change is done in-place!

Parameters
[in]mapPointer to a the internal variable containing the map (it will be changed in place).
[in]xDimIndsThe size of the map x dimension in indices.
[in]yDimIndsThe size of the map y dimension in indices.
[in]zDimIndsThe size of the map z dimension in indices.
[in]weightsA pointer to 1D array of the weights values.
[in]xDimIndsWghThe size of the weights x dimension in indices.
[in]yDimIndsWghThe size of the weights y dimension in indices.
[in]zDimIndsWghThe size of the weights z dimension in indices.
[in]verboseHow much std::out output would you like?
[in]messageShiftAre we in a subprocess, so that the log should be shifted for this function call? If so, by how much?

Definition at line 638 of file ProSHADE_io.cpp.

639 {
640  //================================================ Initialise local variables
641  proshade_double* weightsFinal;
642  size_t origVolume = xDimInds * yDimInds * zDimInds;
643  size_t newVolume = xDimIndsWgh * yDimIndsWgh * zDimIndsWgh;
644 
645  //================================================ If weights have different number of indices than map, then re-sample weights in supplied space
646  if ( ( xDimIndsWgh != xDimInds ) || ( yDimIndsWgh != yDimInds ) || ( zDimIndsWgh != zDimInds ) )
647  {
648  //============================================ Initialise variables
649  fftw_complex* origCoeffs = new fftw_complex [newVolume ];
650  fftw_complex* origCoeffsHKL = new fftw_complex [newVolume ];
651  fftw_complex* modifCoeffs = new fftw_complex [origVolume];
652  fftw_complex* modifCoeffsHKL = new fftw_complex [origVolume];
653  fftw_complex* inMap = new fftw_complex [newVolume ];
654  fftw_complex* outMap = new fftw_complex [origVolume];
655 
656  //============================================ Check memory allocation
657  ProSHADE_internal_misc::checkMemoryAllocation ( inMap, __FILE__, __LINE__, __func__ );
658  ProSHADE_internal_misc::checkMemoryAllocation ( outMap, __FILE__, __LINE__, __func__ );
659  ProSHADE_internal_misc::checkMemoryAllocation ( origCoeffs, __FILE__, __LINE__, __func__ );
660  ProSHADE_internal_misc::checkMemoryAllocation ( modifCoeffs, __FILE__, __LINE__, __func__ );
661  ProSHADE_internal_misc::checkMemoryAllocation ( origCoeffsHKL, __FILE__, __LINE__, __func__ );
662  ProSHADE_internal_misc::checkMemoryAllocation ( modifCoeffsHKL, __FILE__, __LINE__, __func__ );
663 
664  //============================================ Set array to zeroes
665  for ( size_t iter = 0; iter < origVolume; iter++ ) { modifCoeffsHKL[iter][0] = 0.0; modifCoeffsHKL[iter][1] = 0.0; }
666 
667  //============================================ Copy weights to Fourier input array
668  for ( size_t iter = 0; iter < newVolume; iter++ ) { inMap[iter][0] = weights[iter]; inMap[iter][1] = 0.0; }
669 
670  //============================================ Prepare Fourier transform plans
671  fftw_plan planForwardFourier = fftw_plan_dft_3d ( static_cast< int > ( xDimIndsWgh ), static_cast< int > ( yDimIndsWgh ), static_cast< int > ( zDimIndsWgh ), inMap, origCoeffs, FFTW_FORWARD, FFTW_ESTIMATE );
672  fftw_plan inverseFoourier = fftw_plan_dft_3d ( static_cast< int > ( xDimInds ), static_cast< int > ( yDimInds ), static_cast< int > ( zDimInds ), modifCoeffs, outMap, FFTW_BACKWARD, FFTW_ESTIMATE );
673 
674  //============================================ Compute pre and post changes
675  proshade_signed xPre, yPre, zPre;
676  xPre = std::abs ( ( static_cast< proshade_signed > ( xDimInds ) - static_cast< proshade_signed > ( xDimIndsWgh ) ) / 2 );
677  yPre = std::abs ( ( static_cast< proshade_signed > ( yDimInds ) - static_cast< proshade_signed > ( yDimIndsWgh ) ) / 2 );
678  zPre = std::abs ( ( static_cast< proshade_signed > ( zDimInds ) - static_cast< proshade_signed > ( zDimIndsWgh ) ) / 2 );
679 
680  if ( ( ( static_cast< proshade_signed > ( xDimInds ) - static_cast< proshade_signed > ( xDimIndsWgh ) ) % 2 ) == 1 ) { xPre -= 1; }
681  if ( ( ( static_cast< proshade_signed > ( yDimInds ) - static_cast< proshade_signed > ( yDimIndsWgh ) ) % 2 ) == 1 ) { yPre -= 1; }
682  if ( ( ( static_cast< proshade_signed > ( zDimInds ) - static_cast< proshade_signed > ( zDimIndsWgh ) ) % 2 ) == 1 ) { zPre -= 1; }
683 
684  //============================================ Run forward Fourier
685  fftw_execute ( planForwardFourier );
686 
687  //============================================ Initialise local variables
688  proshade_signed maskMapIndex = 0;
689  proshade_signed densMapIndex = 0;
690  proshade_signed xMaskPos, yMaskPos, zMaskPos, xDensPos, yDensPos, zDensPos;
691  proshade_signed maskH, maskK, maskL;
692 
693  //============================================ Convert weights to HKL for re-boxing
694  for ( proshade_signed xIt = 0; xIt < static_cast< proshade_signed > ( xDimIndsWgh ); xIt++ )
695  {
696  for ( proshade_signed yIt = 0; yIt < static_cast< proshade_signed > ( yDimIndsWgh ); yIt++ )
697  {
698  for ( proshade_signed zIt = 0; zIt < static_cast< proshade_signed > ( zDimIndsWgh ); zIt++ )
699  {
700  //================================ Convert to HKL
701  maskH = xIt + static_cast< proshade_signed > ( (xDimIndsWgh+1) / 2 ); if ( maskH >= static_cast< proshade_signed > ( xDimIndsWgh ) ) { maskH -= xDimIndsWgh; }
702  maskK = yIt + static_cast< proshade_signed > ( (yDimIndsWgh+1) / 2 ); if ( maskK >= static_cast< proshade_signed > ( yDimIndsWgh ) ) { maskK -= yDimIndsWgh; }
703  maskL = zIt + static_cast< proshade_signed > ( (zDimIndsWgh+1) / 2 ); if ( maskL >= static_cast< proshade_signed > ( zDimIndsWgh ) ) { maskL -= zDimIndsWgh; }
704 
705  //================================ Find the positions
706  maskMapIndex = zIt + static_cast< proshade_signed > ( zDimIndsWgh ) * ( yIt + static_cast< proshade_signed > ( yDimIndsWgh ) * xIt );
707  densMapIndex = maskL + static_cast< proshade_signed > ( zDimIndsWgh ) * ( maskK + static_cast< proshade_signed > ( yDimIndsWgh ) * maskH );
708 
709  //================================ Save the values
710  origCoeffsHKL[densMapIndex][0] = origCoeffs[maskMapIndex][0];
711  origCoeffsHKL[densMapIndex][1] = origCoeffs[maskMapIndex][1];
712  }
713  }
714  }
715 
716  //============================================ Rebox
717  for ( proshade_signed xIt = 0; xIt < static_cast< proshade_signed > ( xDimInds ); xIt++ )
718  {
719  for ( proshade_signed yIt = 0; yIt < static_cast< proshade_signed > ( yDimInds ); yIt++ )
720  {
721  for ( proshade_signed zIt = 0; zIt < static_cast< proshade_signed > ( zDimInds ); zIt++ )
722  {
723  //================================ Deal with X
724  if ( xDimIndsWgh >= xDimInds ) { xMaskPos = xIt + xPre; }
725  else { xMaskPos = xIt - xPre; }
726  xDensPos = xIt;
727 
728  //================================ Deal with Y
729  if ( yDimIndsWgh >= yDimInds ) { yMaskPos = yIt + yPre; }
730  else { yMaskPos = yIt - yPre; }
731  yDensPos = yIt;
732 
733  //================================ Deal with Z
734  if ( zDimIndsWgh >= zDimInds ) { zMaskPos = zIt + zPre; }
735  else { zMaskPos = zIt - zPre; }
736  zDensPos = zIt;
737 
738  //================================ Skip if weights value not available (because the modifCoeffsHKL array is zeroed, we do not need to do anything here)
739  if ( ( xMaskPos < 0 ) || ( xMaskPos >= static_cast< proshade_signed > ( xDimIndsWgh ) ) ) { continue; }
740  if ( ( yMaskPos < 0 ) || ( yMaskPos >= static_cast< proshade_signed > ( yDimIndsWgh ) ) ) { continue; }
741  if ( ( zMaskPos < 0 ) || ( zMaskPos >= static_cast< proshade_signed > ( zDimIndsWgh ) ) ) { continue; }
742 
743  //================================ Find the positions
744  maskMapIndex = zMaskPos + static_cast< proshade_signed > ( zDimIndsWgh ) * ( yMaskPos + static_cast< proshade_signed > ( yDimIndsWgh ) * xMaskPos );
745  densMapIndex = zDensPos + static_cast< proshade_signed > ( zDimInds ) * ( yDensPos + static_cast< proshade_signed > ( yDimInds ) * xDensPos );
746 
747  //================================ Copy values
748  modifCoeffsHKL[densMapIndex][0] = origCoeffsHKL[maskMapIndex][0];
749  modifCoeffsHKL[densMapIndex][1] = origCoeffsHKL[maskMapIndex][1];
750  }
751  }
752  }
753 
754  //============================================ Convert weights back to FFTW order
755  for ( proshade_signed xIt = 0; xIt < static_cast< proshade_signed > ( xDimInds ); xIt++ )
756  {
757  for ( proshade_signed yIt = 0; yIt < static_cast< proshade_signed > ( yDimInds ); yIt++ )
758  {
759  for ( proshade_signed zIt = 0; zIt < static_cast< proshade_signed > ( zDimInds ); zIt++ )
760  {
761  //================================ Convert to HKL
762  maskH = xIt + static_cast< proshade_signed > ( xDimInds / 2 ); if ( maskH >= static_cast< proshade_signed > ( xDimInds ) ) { maskH -= xDimInds; }
763  maskK = yIt + static_cast< proshade_signed > ( yDimInds / 2 ); if ( maskK >= static_cast< proshade_signed > ( yDimInds ) ) { maskK -= yDimInds; }
764  maskL = zIt + static_cast< proshade_signed > ( zDimInds / 2 ); if ( maskL >= static_cast< proshade_signed > ( zDimInds ) ) { maskL -= zDimInds; }
765 
766  //================================ Find the positions
767  maskMapIndex = zIt + static_cast< proshade_signed > ( zDimInds ) * ( yIt + static_cast< proshade_signed > ( yDimInds ) * xIt );
768  densMapIndex = maskL + static_cast< proshade_signed > ( zDimInds ) * ( maskK + static_cast< proshade_signed > ( yDimInds ) * maskH );
769 
770  //================================ Save the values
771  modifCoeffs[densMapIndex][0] = modifCoeffsHKL[maskMapIndex][0];
772  modifCoeffs[densMapIndex][1] = modifCoeffsHKL[maskMapIndex][1];
773  }
774  }
775  }
776 
777  //============================================ Run inverse Fourier on the modified coefficients
778  fftw_execute ( inverseFoourier );
779 
780  //============================================ Delete old weights and allocate memory for the new, re-sampled weights
781  weightsFinal = new proshade_double [origVolume];
782  ProSHADE_internal_misc::checkMemoryAllocation ( weightsFinal, __FILE__, __LINE__, __func__ );
783 
784  //============================================ Copy results into a new, properly sampled weights
785  for ( size_t iter = 0; iter < origVolume; iter++ ) { weightsFinal[iter] = outMap[iter][0]; }
786 
787  //============================================ Release remaining memory
788  fftw_destroy_plan ( planForwardFourier );
789  fftw_destroy_plan ( inverseFoourier );
790  delete[] origCoeffs;
791  delete[] modifCoeffs;
792  delete[] origCoeffsHKL;
793  delete[] modifCoeffsHKL;
794  delete[] inMap;
795  delete[] outMap;
796  }
797  else
798  {
799  weightsFinal = new proshade_double [origVolume];
800  ProSHADE_internal_misc::checkMemoryAllocation ( weightsFinal, __FILE__, __LINE__, __func__ );
801  for ( size_t iter = 0; iter < origVolume; iter++ ) { weightsFinal[iter] = weights[iter]; }
802  }
803 
804  //================================================ Allocate memory for map Fourier transform
805  fftw_complex* inMap = new fftw_complex [origVolume];
806  fftw_complex* outMap = new fftw_complex [origVolume];
807  ProSHADE_internal_misc::checkMemoryAllocation ( inMap, __FILE__, __LINE__, __func__ );
808  ProSHADE_internal_misc::checkMemoryAllocation ( outMap, __FILE__, __LINE__, __func__ );
809  fftw_plan planForwardFourier = fftw_plan_dft_3d ( static_cast< int > ( xDimInds ), static_cast< int > ( yDimInds ), static_cast< int > ( zDimInds ), inMap, outMap, FFTW_FORWARD, FFTW_ESTIMATE );
810  fftw_plan inverseFoourier = fftw_plan_dft_3d ( static_cast< int > ( xDimInds ), static_cast< int > ( yDimInds ), static_cast< int > ( zDimInds ), outMap, inMap, FFTW_BACKWARD, FFTW_ESTIMATE );
811 
812  //================================================ Set data
813  for ( size_t iter = 0; iter < static_cast< size_t > ( origVolume ); iter++ ) { inMap[iter][0] = map[iter]; inMap[iter][1] = 0.0; }
814 
815  //================================================ Convert map to Fourier space
816  fftw_execute ( planForwardFourier );
817 
818  //================================================ Apply the weights to the map in Fourier space
819  proshade_double normFactor = static_cast<proshade_double> ( origVolume );
820  for ( size_t iter = 0; iter < static_cast< size_t > ( origVolume ); iter++ ) { outMap[iter][0] *= weightsFinal[iter] / normFactor; outMap[iter][1] *= weightsFinal[iter] / normFactor; }
821 
822  //================================================ Convert weighted map from Fourier space
823  fftw_execute ( inverseFoourier );
824 
825  //================================================ Copy results to map
826  for ( size_t iter = 0; iter < static_cast< size_t > ( xDimInds * yDimInds * zDimInds ); iter++ ) { map[iter] = inMap[iter][0]; }
827 
828  //================================================ Release memory
829  delete[] weightsFinal;
830  delete[] inMap;
831  delete[] outMap;
832  fftw_destroy_plan ( planForwardFourier );
833  fftw_destroy_plan ( inverseFoourier );
834 
835  //================================================ Report progress
836  ProSHADE_internal_messages::printProgressMessage ( verbose, 3, "Mask read in and applied successfully.", messageShift );
837 
838  //================================================ Done
839  return ;
840 
841 }

◆ figureDataType()

ProSHADE_internal_io::InputType ProSHADE_internal_io::figureDataType ( std::string  fName)

Function determining input data type.

This function determines the type of the input structure. The possible outputs are MAP for MRC map files, PDB for mmCIF or PDB formatted data, or UNKNOWN if gemmi fail to read the file as co-ordinates as well as map.

Parameters
[in]fNameThe file name of the file for which the type should be determined.
[out]XProSHADE InputType variable with values UNKNOWN, MAP or PDB depending on the type of the input file.

Definition at line 923 of file ProSHADE_io.cpp.

924 {
925  //================================================ Try readin as PDB
926  if ( isFilePDB ( fName ) )
927  {
928  return ( PDB );
929  }
930 
931  //================================================ If not, try readin as MAP
932  if ( isFileMAP ( fName ) )
933  {
934  return ( MAP );
935  }
936 
937  //================================================ No luck? UNKNOWN it is ...
938  return ( UNKNOWN );
939 
940  //================================================ Done
941 
942 }

◆ isFileMAP()

bool ProSHADE_internal_io::isFileMAP ( std::string  fName)

Function determining if the input data type is MAP.

This function checks if the input file is a MAP file and can be read by the CMAP library.

Parameters
[in]fNameThe file name of the file for which the type should be determined.
[out]XBool value true if the file is a MAP file readable by CMAP and false otherwise.

Definition at line 60 of file ProSHADE_io.cpp.

61 {
62  gemmi::Ccp4<float> map;
63  try
64  {
65  map.read_ccp4 ( gemmi::MaybeGzipped (fName.c_str() ) );
66  }
67  catch ( std::runtime_error& e )
68  {
69  //============================================ Supress MSVC C4101 Unferenced variable warning / clang and gcc -Wunused-exception-parameter
70  (void)e;
71 
72  //============================================ Failed to read the map
73  return ( false );
74  }
75 
76  //================================================ Done
77  return ( true );
78 
79 }

◆ isFilePDB()

bool ProSHADE_internal_io::isFilePDB ( std::string  fName)

Function determining if the input data type is PDB.

This function checks if the input file is a PDB file and can be read by the gemmi library.

Parameters
[in]fNameThe file name of the file for which the type should be determined.
[out]XBool value true if the file is a PDB file readable by gemmi and false otherwise.

Definition at line 32 of file ProSHADE_io.cpp.

33 {
34  //================================================ Try reading the file using Gemmi
35  try
36  {
37  gemmi::Structure structure = gemmi::read_structure ( gemmi::MaybeGzipped ( fName ) );
38  }
39  catch ( std::runtime_error& e )
40  {
41  //============================================ Supress MSVC C4101 Unferenced variable warning / clang and gcc -Wunused-exception-parameter
42  (void)e;
43 
44  //============================================ Read failed. Done
45  return ( false );
46  }
47 
48  //================================================ Read successfull. Done
49  return ( true );
50 
51 }

◆ readInMapData() [1/2]

void ProSHADE_internal_io::readInMapData ( gemmi::Ccp4< float > *  gemmiMap,
proshade_double *&  map,
proshade_unsign  xDimInds,
proshade_unsign  yDimInds,
proshade_unsign  zDimInds,
proshade_unsign  xAxOrder,
proshade_unsign  yAxOrder,
proshade_unsign  zAxOrder 
)

This function converts the gemmi Ccp4 object data to ProSHADE internal map representation.

This function firstly allocates the required memory for the ProSHADE internal map representation variable according to the grid size. Then, it iterates over the axes in such a way, so that the resulting ProSHADE variable would have XYZ axis order independently on the axis order of the Ccp4 gemmi object. This should not be necessary as the gemmi setup function should have been called by now, but one never knows.

Parameters
[in]gemmiMapPointer to a gemmi Ccp4 object containing the read in MAP file information.
[in]mapPointer reference to a variable to save the map data.
[in]xDimIndsThe size of x dimension in indices.
[in]yDimIndsThe size of y dimension in indices.
[in]zDimIndsThe size of z dimension in indices.
[in]xAxOrderThe order of the x-axis.
[in]yAxOrderThe order of the y-axis.
[in]zAxOrderThe order of the z-axis.

Definition at line 178 of file ProSHADE_io.cpp.

179 {
180  //================================================ Allocate internal variables
181  proshade_unsign *axOrdArr = new proshade_unsign[3];
182  proshade_unsign *axDimArr = new proshade_unsign[3];
183  proshade_unsign arrPos = 0;
184 
185  //================================================ Check memory allocation and fill in values
186  ProSHADE_internal_misc::checkMemoryAllocation ( axOrdArr, __FILE__, __LINE__, __func__ );
187  ProSHADE_internal_misc::checkMemoryAllocation ( axDimArr, __FILE__, __LINE__, __func__ );
188  axDimArr[0] = xDimInds;
189  axDimArr[1] = yDimInds;
190  axDimArr[2] = zDimInds;
191 
192  //================================================ Allocate the ProSHADE internal map variable memory
193  map = new proshade_double [xDimInds * yDimInds * zDimInds];
194  ProSHADE_internal_misc::checkMemoryAllocation ( map, __FILE__, __LINE__, __func__ );
195 
196  //================================================ Copy read in data to internal map variable
197  for ( axOrdArr[0] = 0; axOrdArr[0] < axDimArr[xAxOrder-1]; axOrdArr[0]++ )
198  {
199  for ( axOrdArr[1] = 0; axOrdArr[1] < axDimArr[yAxOrder-1]; axOrdArr[1]++ )
200  {
201  for ( axOrdArr[2] = 0; axOrdArr[2] < axDimArr[zAxOrder-1]; axOrdArr[2]++ )
202  {
203  arrPos = axOrdArr[2] + axDimArr[zAxOrder-1] * ( axOrdArr[1] + axDimArr[yAxOrder-1] * axOrdArr[0] );
204  map[arrPos] = static_cast< proshade_double > ( gemmiMap->grid.get_value_q( static_cast< int > ( axOrdArr[xAxOrder-1] ),
205  static_cast< int > ( axOrdArr[yAxOrder-1] ),
206  static_cast< int > ( axOrdArr[zAxOrder-1] ) ) );
207  }
208  }
209  }
210 
211  //================================================ Release internal variables memory
212  delete[] axDimArr;
213  delete[] axOrdArr;
214 
215  //================================================ Done
216  return ;
217 
218 }

◆ readInMapData() [2/2]

void ProSHADE_internal_io::readInMapData ( gemmi::Ccp4< int8_t > *  gemmiMap,
proshade_double *&  map,
proshade_unsign  xDimInds,
proshade_unsign  yDimInds,
proshade_unsign  zDimInds,
proshade_unsign  xAxOrder,
proshade_unsign  yAxOrder,
proshade_unsign  zAxOrder 
)

This function converts the gemmi Ccp4 object data to ProSHADE mask representation.

This function firstly allocates the required memory for the ProSHADE mask representation variable according to the grid size. Then, it iterates over the axes in such a way, so that the resulting ProSHADE variable would have XYZ axis order independently on the axis order of the Ccp4 gemmi object. This should not be necessary as the gemmi setup function should have been called by now, but one never knows.

Parameters
[in]gemmiMapPointer to a gemmi Ccp4 object containing the read in mask file information.
[in]mapPointer reference to a variable to save the map data.
[in]xDimIndsThe size of x dimension in indices.
[in]yDimIndsThe size of y dimension in indices.
[in]zDimIndsThe size of z dimension in indices.
[in]xAxOrderThe order of the x-axis.
[in]yAxOrderThe order of the y-axis.
[in]zAxOrderThe order of the z-axis.

Definition at line 234 of file ProSHADE_io.cpp.

235 {
236  //================================================ Allocate internal variables
237  proshade_unsign *axOrdArr = new proshade_unsign[3];
238  proshade_unsign *axDimArr = new proshade_unsign[3];
239  proshade_unsign arrPos = 0;
240 
241  //================================================ Check memory allocation and fill in values
242  ProSHADE_internal_misc::checkMemoryAllocation ( axOrdArr, __FILE__, __LINE__, __func__ );
243  ProSHADE_internal_misc::checkMemoryAllocation ( axDimArr, __FILE__, __LINE__, __func__ );
244  axDimArr[0] = xDimInds;
245  axDimArr[1] = yDimInds;
246  axDimArr[2] = zDimInds;
247 
248  //================================================ Allocate the ProSHADE internal map variable memory
249  map = new proshade_double [xDimInds * yDimInds * zDimInds];
250  ProSHADE_internal_misc::checkMemoryAllocation ( map, __FILE__, __LINE__, __func__ );
251 
252  //================================================ Copy read in data to internal map variable
253  for ( axOrdArr[0] = 0; axOrdArr[0] < axDimArr[xAxOrder-1]; axOrdArr[0]++ )
254  {
255  for ( axOrdArr[1] = 0; axOrdArr[1] < axDimArr[yAxOrder-1]; axOrdArr[1]++ )
256  {
257  for ( axOrdArr[2] = 0; axOrdArr[2] < axDimArr[zAxOrder-1]; axOrdArr[2]++ )
258  {
259  arrPos = axOrdArr[2] + axDimArr[zAxOrder-1] * ( axOrdArr[1] + axDimArr[yAxOrder-1] * axOrdArr[0] );
260  map[arrPos] = static_cast< proshade_double > ( gemmiMap->grid.get_value_q( static_cast< int > ( axOrdArr[xAxOrder-1] ),
261  static_cast< int > ( axOrdArr[yAxOrder-1] ),
262  static_cast< int > ( axOrdArr[zAxOrder-1] ) ) );
263  }
264  }
265  }
266 
267  //================================================ Release internal variables memory
268  delete[] axDimArr;
269  delete[] axOrdArr;
270 
271  //================================================ Done
272  return ;
273 
274 }

◆ readInMapHeader()

void ProSHADE_internal_io::readInMapHeader ( gemmi::Ccp4< float > *  map,
proshade_unsign *  xDimInds,
proshade_unsign *  yDimInds,
proshade_unsign *  zDimInds,
proshade_single *  xDim,
proshade_single *  yDim,
proshade_single *  zDim,
proshade_single *  aAng,
proshade_single *  bAng,
proshade_single *  cAng,
proshade_signed *  xFrom,
proshade_signed *  yFrom,
proshade_signed *  zFrom,
proshade_signed *  xAxOrigin,
proshade_signed *  yAxOrigin,
proshade_signed *  zAxOrigin,
proshade_unsign *  xAxOrder,
proshade_unsign *  yAxOrder,
proshade_unsign *  zAxOrder,
proshade_unsign *  xGridInds,
proshade_unsign *  yGridInds,
proshade_unsign *  zGridInds 
)

This function parses the CCP4 MAP file header as read in by gemmi.

This function uses the gemmi Ccp4 object, which contains all the information read in from a MAP file (including the header), to parse out the ProSHADE required information from the header and saving it to the supplied variables.

Parameters
[in]mapA gemmi Ccp4 objecct containing all the data read in from a MAP file.
[in]xDimIndsAddress to a variable to save the x-axis size in indices.
[in]yDimIndsAddress to a variable to save the y-axis size in indices.
[in]zDimIndsAddress to a variable to save the z-axis size in indices.
[in]xDimAddress to a variable to save the x dimension size in angstroms.
[in]yDimAddress to a variable to save the y dimension size in angstroms.
[in]zDimAddress to a variable to save the z dimension size in angstroms.
[in]aAngAddress to a variable to save the a angle in degrees.
[in]bAngAddress to a variable to save the b angle in degrees.
[in]cAngAddress to a variable to save the c angle in degrees.
[in]xFromAddress to a variable to save the starting index along the x-axis.
[in]yFromAddress to a variable to save the starting index along the y-axis.
[in]zFromAddress to a variable to save the starting index along the z-axis.
[in]xAxOriginAddress to a variable to save the map origin positon along the x-axis.
[in]yAxOriginAddress to a variable to save the map origin positon along the y-axis.
[in]zAxOriginAddress to a variable to save the map origin positon along the z-axis.
[in]xAxOrderAddress to a variable to save the order of x axis.
[in]yAxOrderAddress to a variable to save the order of y axis.
[in]zAxOrderAddress to a variable to save the order of z axis.
[in]xGridIndsAddress to a variable to save the grid indices count along the x-axis.
[in]yGridIndsAddress to a variable to save the grid indices count along the y-axis.
[in]zGridIndsAddress to a variable to save the grid indices count along the z-axis.

Definition at line 109 of file ProSHADE_io.cpp.

110 {
111  //================================================ Read in the map file header
112  *xDimInds = static_cast<proshade_unsign> ( map->header_i32 ( 1 ) );
113  *yDimInds = static_cast<proshade_unsign> ( map->header_i32 ( 2 ) );
114  *zDimInds = static_cast<proshade_unsign> ( map->header_i32 ( 3 ) );
115 
116  *xFrom = static_cast<proshade_signed> ( map->header_i32 ( 5 ) );
117  *yFrom = static_cast<proshade_signed> ( map->header_i32 ( 6 ) );
118  *zFrom = static_cast<proshade_signed> ( map->header_i32 ( 7 ) );
119 
120  *xDim = static_cast<proshade_single> ( map->header_float ( 11 ) );
121  *yDim = static_cast<proshade_single> ( map->header_float ( 12 ) );
122  *zDim = static_cast<proshade_single> ( map->header_float ( 13 ) );
123 
124  *aAng = static_cast<proshade_single> ( map->header_float ( 14 ) );
125  *bAng = static_cast<proshade_single> ( map->header_float ( 15 ) );
126  *cAng = static_cast<proshade_single> ( map->header_float ( 16 ) );
127 
128  *xAxOrigin = static_cast<proshade_signed> ( map->header_i32 ( 50 ) ) + (*xFrom);
129  *yAxOrigin = static_cast<proshade_signed> ( map->header_i32 ( 51 ) ) + (*yFrom);
130  *zAxOrigin = static_cast<proshade_signed> ( map->header_i32 ( 52 ) ) + (*zFrom);
131 
132  *xAxOrder = static_cast<proshade_unsign> ( map->header_i32 ( 17 ) );
133  *yAxOrder = static_cast<proshade_unsign> ( map->header_i32 ( 18 ) );
134  *zAxOrder = static_cast<proshade_unsign> ( map->header_i32 ( 19 ) );
135 
136  *xGridInds = static_cast<proshade_unsign> ( map->header_i32 ( 8 ) );
137  *yGridInds = static_cast<proshade_unsign> ( map->header_i32 ( 9 ) );
138  *zGridInds = static_cast<proshade_unsign> ( map->header_i32 ( 10 ) );
139 
140  //================================================ Deal with sampling being different from cell size
141  if ( *xGridInds != *xDimInds )
142  {
143  *xDim = *xDim * ( static_cast<proshade_single> ( *xDimInds ) / static_cast<proshade_single> ( *xGridInds ) );
144  *xGridInds = *xDimInds;
145  }
146 
147  if ( *yGridInds != *yDimInds )
148  {
149  *yDim = *yDim * ( static_cast<proshade_single> ( *yDimInds ) / static_cast<proshade_single> ( *yGridInds ) );
150  *yGridInds = *yDimInds;
151  }
152 
153  if ( *zGridInds != *zDimInds )
154  {
155  *zDim = *zDim * ( static_cast<proshade_single> ( *zDimInds ) / static_cast<proshade_single> ( *zGridInds ) );
156  *zGridInds = *zDimInds;
157  }
158 
159  //================================================ Done
160  return ;
161 
162 }

◆ writeOutMapHeader()

void ProSHADE_internal_io::writeOutMapHeader ( gemmi::Ccp4< float > *  map,
proshade_unsign  xDimInds,
proshade_unsign  yDimInds,
proshade_unsign  zDimInds,
proshade_single  xDim,
proshade_single  yDim,
proshade_single  zDim,
proshade_single  aAng,
proshade_single  bAng,
proshade_single  cAng,
proshade_signed  xFrom,
proshade_signed  yFrom,
proshade_signed  zFrom,
proshade_signed  xAxOrigin,
proshade_signed  yAxOrigin,
proshade_signed  zAxOrigin,
proshade_unsign  xAxOrder,
proshade_unsign  yAxOrder,
proshade_unsign  zAxOrder,
proshade_unsign  xGridInds,
proshade_unsign  yGridInds,
proshade_unsign  zGridInds,
std::string  title,
int  mode 
)

This function parses the CCP4 MAP file header as read in by gemmi.

This function uses the gemmi Ccp4 object, which contains all the information read in from a MAP file (including the header), to parse out the ProSHADE required information from the header and saving it to the supplied variables.

Parameters
[in]mapA gemmi Ccp4 objecct containing all the data read in from a MAP file.
[in]xDimIndsVariable holding the x-axis size in indices.
[in]yDimIndsVariable holding the y-axis size in indices.
[in]zDimIndsVariable holding the z-axis size in indices.
[in]xDimVariable holding the x dimension size in angstroms.
[in]yDimVariable holding the y dimension size in angstroms.
[in]zDimVariable holding the z dimension size in angstroms.
[in]aAngVariable holding the a angle in degrees.
[in]bAngVariable holding the b angle in degrees.
[in]cAngVariable holding the c angle in degrees.
[in]xFromVariable holding the starting index along the x-axis.
[in]yFromVariable holding the starting index along the y-axis.
[in]zFromVariable holding the starting index along the z-axis.
[in]xAxOriginVariable holding the map origin positon along the x-axis.
[in]yAxOriginVariable holding the map origin positon along the y-axis.
[in]zAxOriginVariable holding the map origin positon along the z-axis.
[in]xAxOrderVariable holding the order of x axis.
[in]yAxOrderVariable holding the order of y axis.
[in]zAxOrderVariable holding the order of z axis.
[in]xGridIndsVariable holding the grid indices count along the x-axis.
[in]yGridIndsVariable holding the grid indices count along the y-axis.
[in]zGridIndsVariable holding the grid indices count along the z-axis.
[in]titleThe title to be written into the MAP file.
[in]modeThe variable type of the data, please leave two (float) unless you require any specific other mode.

Definition at line 873 of file ProSHADE_io.cpp.

874 {
875  //================================================ Fill in the map file header
876  map->set_header_i32 ( 1 , static_cast<int32_t> ( xDimInds ) ); // Number of columns in 3D data array (fast axis)
877  map->set_header_i32 ( 2 , static_cast<int32_t> ( yDimInds ) ); // Number of columns in 3D data array (medium axis)
878  map->set_header_i32 ( 3 , static_cast<int32_t> ( zDimInds ) ); // Number of columns in 3D data array (slow axis)
879  map->set_header_i32 ( 4 , static_cast<int32_t> ( mode ) ); // Map mode
880  map->set_header_i32 ( 5 , static_cast<int32_t> ( xFrom ) ); // Starting index (fast axis)
881  map->set_header_i32 ( 6 , static_cast<int32_t> ( yFrom ) ); // Starting index (medium axis)
882  map->set_header_i32 ( 7 , static_cast<int32_t> ( zFrom ) ); // Starting index (slow axis)
883  map->set_header_i32 ( 8 , static_cast<int32_t> ( xGridInds ) ); // Grid sampling (fast axis)
884  map->set_header_i32 ( 9 , static_cast<int32_t> ( yGridInds ) ); // Grid sampling (medium axis)
885  map->set_header_i32 ( 10, static_cast<int32_t> ( zGridInds ) ); // Grid sampling (slow axis)
886  map->set_header_float ( 11, static_cast<float> ( xDim ) ); // Grid dimension in Angstrom (fast axis)
887  map->set_header_float ( 12, static_cast<float> ( yDim ) ); // Grid dimension in Angstrom (medium axis)
888  map->set_header_float ( 13, static_cast<float> ( zDim ) ); // Grid dimension in Angstrom (slow axis)
889  map->set_header_float ( 14, static_cast<float> ( aAng ) ); // Alpha angle in degrees
890  map->set_header_float ( 15, static_cast<float> ( bAng ) ); // Beta angle in degrees
891  map->set_header_float ( 16, static_cast<float> ( cAng ) ); // Gamma angle in degrees
892  map->set_header_i32 ( 17, static_cast<int32_t> ( xAxOrder ) ); // MAPC
893  map->set_header_i32 ( 18, static_cast<int32_t> ( yAxOrder ) ); // MAPR
894  map->set_header_i32 ( 19, static_cast<int32_t> ( zAxOrder ) ); // MAPS
895  if ( map->grid.spacegroup ) { map->set_header_i32 ( 23, static_cast<int32_t> ( map->grid.spacegroup->ccp4 ) ); } // Space group
896  else { map->set_header_i32 ( 23, static_cast<int32_t> ( 1 ) ); }
897  map->set_header_i32 ( 24, static_cast<int32_t> ( map->grid.spacegroup->operations().order() * 80 ) ); // NSYMBT - size of extended header (which follows main header) in bytes
898  map->set_header_str ( 27, "CCP4" ); // Code for the type of extended header
899  map->set_header_i32 ( 28, static_cast<int32_t> ( 20140 ) ); // Version
900  map->set_header_i32 ( 50, static_cast<int32_t> ( xAxOrigin ) ); // Origin of the map (fast axis)
901  map->set_header_i32 ( 51, static_cast<int32_t> ( yAxOrigin ) ); // Origin of the map (medium axis)
902  map->set_header_i32 ( 52, static_cast<int32_t> ( zAxOrigin ) ); // Origin of the map (slow axis)
903  map->set_header_str ( 53, "MAP" ); // File format
904  if ( gemmi::is_little_endian() ) { map->set_header_i32 ( 54, static_cast<int32_t> ( 0x00004144 ) ); } // Machine stamp encoding byte ordering of data
905  else { map->set_header_i32 ( 54, static_cast<int32_t> ( 0x11110000 ) ); }
906  map->set_header_i32 ( 56, static_cast<int32_t> ( 1 ) ); // Number of labels used
907  std::memset ( reinterpret_cast<void*> ( &(map->ccp4_header.at( 56 )) ), ' ', static_cast< size_t > ( 800 + map->grid.spacegroup->operations().order() * 80 ) ); // 56 is used because the vector is indexed from 0
908  map->set_header_str ( 57, title ); // Title
909 
910  //================================================ Done
911  return ;
912 
913 }

◆ writeRotationTranslationJSON()

void ProSHADE_internal_io::writeRotationTranslationJSON ( proshade_double  trsX1,
proshade_double  trsY1,
proshade_double  trsZ1,
proshade_double  eulA,
proshade_double  eulB,
proshade_double  eulG,
proshade_double  trsX2,
proshade_double  trsY2,
proshade_double  trsZ2,
std::string  fileName 
)

Function for writing out the optimal rotation and translation into a JSON file.

This function takes the initial translation (assuming that the centre of rotation is not at the centre of indices), the rotation (in terms of the Euler angles) and the translation detected by the overlay task and proceeds to write a JSON file containing all of these information in the order in which they should be applied with the exception that depending on around which point the rotation should be done, either the translation to origin or to the map centre should be applied.

This function assumes that the second translation includes the reverse of the first translation already.

Parameters
[in]trsX1The translation required to get the rotation centre to origin along the x-axis.
[in]trsY1The translation required to get the rotation centre to origin along the y-axis.
[in]trsZ1The translation required to get the rotation centre to origin along the z-axis.
[in]eulAThe optimal rotation Euler angle alpha.
[in]eulBThe optimal rotation Euler angle beta.
[in]eulGThe optimal rotation Euler angle gamma.
[in]trsX2The optimal translation along the x-axis + reverse of the trsX1.
[in]trsY2The optimal translation along the y-axis + reverse of the trsY1.
[in]trsZ2The optimal translation along the z-axis + reverse of the trsZ1.
[in]fileNameThe file name of the file for which the information should be written into.

Definition at line 965 of file ProSHADE_io.cpp.

966 {
967  //================================================ Open file for writing
968  std::ofstream jsonFile;
969  jsonFile.open ( fileName );
970 
971  //================================================ Check file opening success
972  if ( !jsonFile.is_open( ) )
973  {
974  throw ProSHADE_exception ( "Failed to open JSON output file.", "E000056", __FILE__, __LINE__, __func__, "Failed to open json file to which the rotation and\n : translation would be written into. Most likely cause is\n : lack of rights to write in the current folder." );
975  }
976 
977  //================================================ Get rotation matrix from Euler angles
978  proshade_double* rotMat = new proshade_double[9];
979  ProSHADE_internal_misc::checkMemoryAllocation ( rotMat, __FILE__, __LINE__, __func__ );
981 
982  //================================================ Write the info
983  jsonFile << "{\n";
984  jsonFile << " \"translationToOrigin\" : [ " << trsX1 << ", " << trsY1 << ", " << trsZ1 << " ], \n";
985 
986  jsonFile << " \"rotationMatrix:\" : [ " << rotMat[0] << ", " << rotMat[1] << ", " << rotMat[2] << ", \n";
987  jsonFile << " " << rotMat[3] << ", " << rotMat[4] << ", " << rotMat[5] << ", \n";
988  jsonFile << " " << rotMat[6] << ", " << rotMat[7] << ", " << rotMat[8] << "], \n";
989 
990  jsonFile << " \"translationFromRotCenToOverlay\" : [ " << trsX2 << ", " << trsY2 << ", " << trsZ2 << " ] \n";
991  jsonFile << "}\n";
992 
993  //================================================ Close file
994  jsonFile.close ( );
995 
996  //================================================ Release memory
997  delete[] rotMat;
998 
999  //================================================ Done
1000  return ;
1001 
1002 }
ProSHADE_internal_io::isFilePDB
bool isFilePDB(std::string fName)
Function determining if the input data type is PDB.
Definition: ProSHADE_io.cpp:32
ProSHADE_exception
This class is the representation of ProSHADE exception.
Definition: ProSHADE_exceptions.hpp:37
ProSHADE_internal_io::applyWeightsFromArray
void applyWeightsFromArray(proshade_double *&map, proshade_unsign xDimInds, proshade_unsign yDimInds, proshade_unsign zDimInds, proshade_double *&weights, proshade_unsign xDimIndsWgh, proshade_unsign yDimIndsWgh, proshade_unsign zDimIndsWgh, proshade_signed verbose, proshade_signed messageShift)
This function applies the weights to the map.
Definition: ProSHADE_io.cpp:638
ProSHADE_internal_maths::getRotationMatrixFromEulerZXZAngles
void getRotationMatrixFromEulerZXZAngles(proshade_double eulerAlpha, proshade_double eulerBeta, proshade_double eulerGamma, proshade_double *matrix)
Function to find the rotation matrix from Euler angles (ZXZ convention).
Definition: ProSHADE_maths.cpp:1020
ProSHADE_internal_io::isFileMAP
bool isFileMAP(std::string fName)
Function determining if the input data type is MAP.
Definition: ProSHADE_io.cpp:60
ProSHADE_internal_io::readInMapData
void readInMapData(gemmi::Ccp4< float > *gemmiMap, proshade_double *&map, proshade_unsign xDimInds, proshade_unsign yDimInds, proshade_unsign zDimInds, proshade_unsign xAxOrder, proshade_unsign yAxOrder, proshade_unsign zAxOrder)
This function converts the gemmi Ccp4 object data to ProSHADE internal map representation.
Definition: ProSHADE_io.cpp:178
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:68
ProSHADE_internal_io::applyMaskFromArray
void applyMaskFromArray(proshade_double *&map, proshade_unsign xDimInds, proshade_unsign yDimInds, proshade_unsign zDimInds, proshade_double *&mask, proshade_unsign xDimIndsMsk, proshade_unsign yDimIndsMsk, proshade_unsign zDimIndsMsk, proshade_signed verbose, proshade_signed messageShift)
This function applies the mask to the map.
Definition: ProSHADE_io.cpp:367
ProSHADE_internal_messages::printProgressMessage
void printProgressMessage(proshade_signed verbose, proshade_signed messageLevel, std::string message, proshade_signed messageShift=0)
General stdout message printing.
Definition: ProSHADE_messages.cpp:71
ProSHADE_internal_io::readInMapHeader
void readInMapHeader(gemmi::Ccp4< float > *map, proshade_unsign *xDimInds, proshade_unsign *yDimInds, proshade_unsign *zDimInds, proshade_single *xDim, proshade_single *yDim, proshade_single *zDim, proshade_single *aAng, proshade_single *bAng, proshade_single *cAng, proshade_signed *xFrom, proshade_signed *yFrom, proshade_signed *zFrom, proshade_signed *xAxOrigin, proshade_signed *yAxOrigin, proshade_signed *zAxOrigin, proshade_unsign *xAxOrder, proshade_unsign *yAxOrder, proshade_unsign *zAxOrder, proshade_unsign *xGridInds, proshade_unsign *yGridInds, proshade_unsign *zGridInds)
This function parses the CCP4 MAP file header as read in by gemmi.
Definition: ProSHADE_io.cpp:109