Xmipp  v3.23.11-Nereus
Public Member Functions | Public Attributes | List of all members

#include <mpi_classify_CLTomo.h>

Inheritance diagram for ProgClassifyCL3D:
Inheritance graph
[legend]
Collaboration diagram for ProgClassifyCL3D:
Collaboration graph
[legend]

Public Member Functions

 ProgClassifyCL3D (int argc, char **argv)
 MPI constructor. More...
 
 ProgClassifyCL3D (const ProgClassifyCL3D &)=delete
 
 ProgClassifyCL3D (const ProgClassifyCL3D &&)=delete
 
 ~ProgClassifyCL3D ()
 Destructor. More...
 
ProgClassifyCL3Doperator= (const ProgClassifyCL3D &)=delete
 
ProgClassifyCL3Doperator= (const ProgClassifyCL3D &&)=delete
 
void readParams ()
 Read. More...
 
void show () const
 Show. More...
 
void defineParams ()
 Usage. More...
 
void produceSideInfo ()
 Produce side info. More...
 
void run ()
 Run. More...
 
- Public Member Functions inherited from XmippProgram
const char * getParam (const char *param, int arg=0)
 
const char * getParam (const char *param, const char *subparam, int arg=0)
 
int getIntParam (const char *param, int arg=0)
 
int getIntParam (const char *param, const char *subparam, int arg=0)
 
double getDoubleParam (const char *param, int arg=0)
 
double getDoubleParam (const char *param, const char *subparam, int arg=0)
 
float getFloatParam (const char *param, int arg=0)
 
float getFloatParam (const char *param, const char *subparam, int arg=0)
 
void getListParam (const char *param, StringVector &list)
 
int getCountParam (const char *param)
 
bool checkParam (const char *param)
 
bool existsParam (const char *param)
 
void addParamsLine (const String &line)
 
void addParamsLine (const char *line)
 
ParamDefgetParamDef (const char *param) const
 
virtual void quit (int exit_code=0) const
 
virtual int tryRun ()
 
void initProgress (size_t total, size_t stepBin=60)
 
void setProgress (size_t value=0)
 
void endProgress ()
 
void processDefaultComment (const char *param, const char *left)
 
void setDefaultComment (const char *param, const char *comment)
 
virtual void initComments ()
 
void setProgramName (const char *name)
 
void addUsageLine (const char *line, bool verbatim=false)
 
void clearUsage ()
 
void addExampleLine (const char *example, bool verbatim=true)
 
void addSeeAlsoLine (const char *seeAlso)
 
void addKeywords (const char *keywords)
 
const char * name () const
 
virtual void usage (int verb=0) const
 
virtual void usage (const String &param, int verb=2)
 
int version () const
 
virtual void read (int argc, const char **argv, bool reportErrors=true)
 
virtual void read (int argc, char **argv, bool reportErrors=true)
 
void read (const String &argumentsLine)
 
 XmippProgram ()
 
 XmippProgram (int argc, const char **argv)
 
virtual ~XmippProgram ()
 

Public Attributes

FileName fnSel
 Input selfile with the images to quantify. More...
 
FileName fnCodes0
 Input selfile with initial codes. More...
 
FileName fnOut
 Output rootname. More...
 
FileName fnSym
 Symmetry file or code. More...
 
int Niter
 Number of iterations. More...
 
int Ncodes0
 Initial number of code vectors. More...
 
int Ncodes
 Final number of code vectors. More...
 
int Nneighbours
 Number of neighbours. More...
 
double PminSize
 Minimum size of a node. More...
 
double maxFreq
 Maximum frequency for the alignment. More...
 
double sparsity
 Sparsity factor (0<f<1; 1=drop all coefficients, 0=do not drop any coefficient) More...
 
double DWTsparsity
 DWT Sparsity factor (0<f<1; 1=drop all coefficients, 0=do not drop any coefficient) More...
 
bool classifyAllImages
 Clasify all images. More...
 
bool randomizeStartingOrientation
 Use this option to avoid aligning at the beginning all the missing wedges. More...
 
double maxShiftZ
 Maximum shift Z. More...
 
double maxShiftY
 Maximum shift Y. More...
 
double maxShiftX
 Maximum shift X. More...
 
double maxRot
 Maximum rot. More...
 
double maxTilt
 Maximum tilt. More...
 
double maxPsi
 Maximum psi. More...
 
Mask mask
 Mask. More...
 
bool dontAlign
 Don't align. More...
 
bool generateAlignedVolumes
 Generate aligned subvolumes. More...
 
SymList SL
 
MetaDataDb SF
 
std::vector< size_t > objId
 
CL3D vq
 
MpiNodenode
 
GaussianInterpolator gaussianInterpolator
 
size_t Zdim
 
size_t Ydim
 
size_t Xdim
 
PyObject * frmFunc
 Pointer to the Python FRM alignment function. More...
 
PyObject * wedgeClass
 Pointer to the Python GeneralWedge class. More...
 
double maxShift
 Max shift. More...
 
MultidimArray< unsigned char > maxFreqMask
 MaxFreq mask. More...
 
- Public Attributes inherited from XmippProgram
bool doRun
 
bool runWithoutArgs
 
int verbose
 Verbosity level. More...
 
int debug
 

Additional Inherited Members

- Protected Member Functions inherited from XmippProgram
void defineCommons ()
 
- Protected Attributes inherited from XmippProgram
int errorCode
 
ProgramDefprogDef
 Program definition and arguments parser. More...
 
std::map< String, CommentListdefaultComments
 
int argc
 Original command line arguments. More...
 
const char ** argv
 

Detailed Description

CL3D parameters.

Definition at line 220 of file mpi_classify_CLTomo.h.

Constructor & Destructor Documentation

◆ ProgClassifyCL3D() [1/3]

ProgClassifyCL3D::ProgClassifyCL3D ( int  argc,
char **  argv 
)

MPI constructor.

Definition at line 1300 of file mpi_classify_CLTomo_prog.cpp.

1301 {
1302  node = new MpiNode(argc, argv);
1303  if (!node->isMaster())
1304  verbose = 0;
1306 }
int allowed_data_types
Definition: mask.h:495
int argc
Original command line arguments.
Definition: xmipp_program.h:86
const char ** argv
Definition: xmipp_program.h:87
int verbose
Verbosity level.
#define INT_MASK
Definition: mask.h:385
bool isMaster() const
Definition: xmipp_mpi.cpp:166

◆ ProgClassifyCL3D() [2/3]

ProgClassifyCL3D::ProgClassifyCL3D ( const ProgClassifyCL3D )
delete

◆ ProgClassifyCL3D() [3/3]

ProgClassifyCL3D::ProgClassifyCL3D ( const ProgClassifyCL3D &&  )
delete

◆ ~ProgClassifyCL3D()

ProgClassifyCL3D::~ProgClassifyCL3D ( )

Destructor.

Definition at line 1309 of file mpi_classify_CLTomo_prog.cpp.

1310 {
1311  delete node;
1312 }

Member Function Documentation

◆ defineParams()

void ProgClassifyCL3D::defineParams ( )
virtual

Usage.

Reimplemented from XmippProgram.

Definition at line 1377 of file mpi_classify_CLTomo_prog.cpp.

1378 {
1379  addUsageLine("Divide a selfile of volumes into the desired number of classes. ");
1380  addUsageLine("+Vector quantization with correntropy and a probabilistic criterion is used for creating the subdivisions.");
1381  addUsageLine("+Correlation and the standard maximum correlation criterion can also be used and normally produce good results.");
1382  addUsageLine("+Correntropy and the probabilistic clustering criterion are recommended for images with very low SNR or cases in which the correlation have clear difficulties to converge.");
1383  addUsageLine("+");
1384  addUsageLine("+The algorithm is fully described in [[http://www.ncbi.nlm.nih.gov/pubmed/20362059][this article]].");
1385  addUsageLine("+");
1386  addUsageLine("+An interesting convergence criterion is the number of images changing classes between iterations. If a low percentage of the image change class, then the clustering is rather stable and clear.");
1387  addUsageLine("+If many images change class, it is likely that there is not enough SNR to determine so many classes. It is recommended to reduce the number of classes");
1388  addParamsLine(" -i <selfile> : Selfile with the input images");
1389  addParamsLine(" [--oroot <root=class>] : Output rootname, by default, class");
1390  addParamsLine(" [--iter <N=20>] : Number of iterations");
1391  addParamsLine(" [--nref0 <N=2>] : Initial number of code vectors");
1392  addParamsLine("or --ref0 <selfile=\"\"> : Selfile with initial code vectors");
1393  addParamsLine(" [--nref <N=16>] : Final number of code vectors");
1394  addParamsLine(" [--neigh+ <N=4>] : Number of neighbour code vectors");
1395  addParamsLine(" : Set -1 for all");
1396  addParamsLine(" [--minsize+ <N=20>] : Percentage minimum node size");
1397  addParamsLine(" [--sparsity+ <f=0.975>] : Percentage of Fourier coefficients to drop");
1398  addParamsLine(" [--DWTsparsity+ <f=0.99>] : Percentage of wavelet coefficients to drop");
1399  addParamsLine(" [--maxShiftX <d=10>] : Maximum allowed shift in X");
1400  addParamsLine(" [--maxShiftY <d=10>] : Maximum allowed shift in Y");
1401  addParamsLine(" [--maxShiftZ <d=10>] : Maximum allowed shift in Z");
1402  addParamsLine(" [--maxRot <d=360>] : Maximum allowed rotational angle in absolute value");
1403  addParamsLine(" [--maxTilt <d=360>] : Maximum allowed tilt angle in absolute value");
1404  addParamsLine(" [--maxPsi <d=360>] : Maximum allowed in-plane angle in absolute value");
1405  addParamsLine(" [--classifyAllImages] : By default, some images may not be classified. Use this option to classify them all.");
1406  addParamsLine(" [--sym <s=c1>] : Symmetry of the classes to be reconstructed");
1407  addParamsLine(" [--maxFreq <w=0.2>] : Maximum frequency to be reconstructed");
1408  addParamsLine(" [--randomizeStartingOrientation] : Use this option to avoid aligning all missing wedges");
1409  addParamsLine(" [--dontAlign] : Do not align volumes, only classify");
1410  addParamsLine(" [--generateAlignedVolumes]: Generate aligned subvolumes at the end");
1411  Mask::defineParams(this,INT_MASK,NULL,NULL,true);
1412  addExampleLine("The MPI program as to be called through a python wrapper that encapsulates some path setting",false);
1413  addExampleLine(" ",false);
1414  addExampleLine("Call the program as xmipp_classify_CLTomo [numberOfMPIProcessors] [arguments as above]",false);
1415  addExampleLine("xmipp_classify_CLTomo 8 -i images.stk --nref 256 --oroot class --iter 10");
1416 }
static void defineParams(XmippProgram *program, int allowed_data_types=ALL_KINDS, const char *prefix=nullptr, const char *comment=nullptr, bool moreOptions=false)
Definition: mask.cpp:1203
void addExampleLine(const char *example, bool verbatim=true)
#define INT_MASK
Definition: mask.h:385
void addUsageLine(const char *line, bool verbatim=false)
void addParamsLine(const String &line)

◆ operator=() [1/2]

ProgClassifyCL3D& ProgClassifyCL3D::operator= ( const ProgClassifyCL3D )
delete

◆ operator=() [2/2]

ProgClassifyCL3D& ProgClassifyCL3D::operator= ( const ProgClassifyCL3D &&  )
delete

◆ produceSideInfo()

void ProgClassifyCL3D::produceSideInfo ( )

Produce side info.

Definition at line 1418 of file mpi_classify_CLTomo_prog.cpp.

1419 {
1420  gaussianInterpolator.initialize(6, 60000, false);
1421 
1422  // Get image dimensions
1423  SF.read(fnSel);
1424  size_t Ndim;
1425  getImageSize(SF, Xdim, Ydim, Zdim, Ndim);
1426 
1427  // MaxFreq mask
1430  V.initZeros(Zdim,Ydim,Xdim);
1431  FourierTransformer transformer;
1432  transformer.FourierTransform(V,VF);
1433  maxFreqMask.initZeros(VF);
1434  Matrix1D<int> idx(3);
1435  Matrix1D<double> freq(3);
1437  {
1438  XX(idx)=j;
1439  YY(idx)=i;
1440  ZZ(idx)=k;
1441  FFT_idx2digfreq(V,idx,freq);
1442  if (freq.module()<=maxFreq)
1443  A3D_ELEM(maxFreqMask,k,i,j)=1;
1444  }
1445 
1446  // Prepare symmetry list
1448 
1449  // Prepare the Task distributor
1450  SF.findObjects(objId);
1451 
1452  // Prepare mask for evaluating the noise outside
1453  MultidimArray<int> sphericalMask;
1454  sphericalMask.resize(Zdim, Ydim, Xdim);
1455  sphericalMask.setXmippOrigin();
1456  BinaryCircularMask(sphericalMask, Xdim / 2, INNER_MASK);
1457 
1460  mMask.setXmippOrigin();
1461  FOR_ALL_ELEMENTS_IN_ARRAY3D(sphericalMask)
1462  if (A3D_ELEM(sphericalMask,k,i,j)==0)
1463  A3D_ELEM(mMask,k,i,j)=0;
1464 
1465  // Read input code vectors if available
1466  std::vector<MultidimArray<double> > codes0;
1467  if (fnCodes0 != "")
1468  {
1469  Image<double> I;
1470  MetaDataDb SFCodes(fnCodes0);
1471 
1472  for (size_t objId : SFCodes.ids())
1473  {
1474  I.readApplyGeo(SFCodes, objId);
1475  I().setXmippOrigin();
1476  codes0.push_back(I());
1477  }
1478  Ncodes0 = codes0.size();
1479  }
1480  vq.initialize(SF, codes0);
1481 }
GaussianInterpolator gaussianInterpolator
void resize(size_t Ndim, size_t Zdim, size_t Ydim, size_t Xdim, bool copy=true)
std::vector< size_t > objId
void findObjects(std::vector< size_t > &objectsOut, const MDQuery &query) const override
void BinaryCircularMask(MultidimArray< int > &mask, double radius, int mode, double x0, double y0, double z0)
Definition: mask.cpp:193
int readSymmetryFile(FileName fn_sym, double accuracy=SYM_ACCURACY)
Definition: symmetries.cpp:33
void getImageSize(const MetaData &md, size_t &Xdim, size_t &Ydim, size_t &Zdim, size_t &Ndim, MDLabel image_label)
void initialize(MetaDataDb &_SF, std::vector< MultidimArray< double > > &_codes0)
Initialize.
FileName fnSel
Input selfile with the images to quantify.
void FFT_idx2digfreq(T &v, const Matrix1D< int > &idx, Matrix1D< double > &freq)
Definition: xmipp_fft.h:80
int readApplyGeo(const FileName &name, const MDRow &row, const ApplyGeoParams &params=DefaultApplyGeoParams)
#define i
ql0001_ & k(htemp+1),(cvec+1),(atemp+1),(bj+1),(bl+1),(bu+1),(x+1),(clamda+1), &iout, infoqp, &zero,(w+1), &lenw,(iw+1), &leniw, &glob_grd.epsmac
#define A3D_ELEM(V, k, i, j)
#define FOR_ALL_ELEMENTS_IN_ARRAY3D(V)
#define XX(v)
Definition: matrix1d.h:85
void FourierTransform(T &v, T1 &V, bool getCopy=true)
Definition: xmipp_fftw.h:166
#define j
#define YY(v)
Definition: matrix1d.h:93
MultidimArray< unsigned char > maxFreqMask
MaxFreq mask.
void generate_mask(bool apply_geo=false)
Definition: mask.cpp:1577
void initialize(double xmax, int N, bool normalize=true)
void read(const FileName &inFile, const std::vector< MDLabel > *desiredLabels=NULL, bool decomposeStack=true) override
void initZeros(const MultidimArray< T1 > &op)
int Ncodes0
Initial number of code vectors.
const MultidimArray< int > & get_binary_mask() const
Definition: mask.h:707
FileName fnCodes0
Input selfile with initial codes.
FileName fnSym
Symmetry file or code.
#define ZZ(v)
Definition: matrix1d.h:101
double maxFreq
Maximum frequency for the alignment.
constexpr int INNER_MASK
Definition: mask.h:47

◆ readParams()

void ProgClassifyCL3D::readParams ( )
virtual

Read.

Reimplemented from XmippProgram.

Definition at line 1315 of file mpi_classify_CLTomo_prog.cpp.

1316 {
1317  fnSel = getParam("-i");
1318  fnOut = getParam("--oroot");
1319  fnCodes0 = getParam("--ref0");
1320  Niter = getIntParam("--iter");
1321  Nneighbours = getIntParam("--neigh");
1322  Ncodes0 = getIntParam("--nref0");
1323  Ncodes = getIntParam("--nref");
1324  PminSize = getDoubleParam("--minsize");
1325  sparsity = getDoubleParam("--sparsity");
1326  DWTsparsity = getDoubleParam("--DWTsparsity");
1327  maxShiftZ = getDoubleParam("--maxShiftZ");
1328  maxShiftY = getDoubleParam("--maxShiftY");
1329  maxShiftX = getDoubleParam("--maxShiftX");
1332  maxFreq=getDoubleParam("--maxFreq");
1333  maxRot = getDoubleParam("--maxRot");
1334  maxTilt = getDoubleParam("--maxTilt");
1335  maxPsi = getDoubleParam("--maxPsi");
1336  classifyAllImages = checkParam("--classifyAllImages");
1337  randomizeStartingOrientation = checkParam("--randomizeStartingOrientation");
1338  fnSym = getParam("--sym");
1339  if (checkParam("--mask"))
1340  mask.readParams(this);
1341  dontAlign = checkParam("--dontAlign");
1342  generateAlignedVolumes = checkParam("--generateAlignedVolumes");
1343 
1344  prmCL3Dprog = this; // FIXME HACK because of the global variable. Solve it properly
1345 }
double maxShift
Max shift.
double sparsity
Sparsity factor (0<f<1; 1=drop all coefficients, 0=do not drop any coefficient)
ProgClassifyCL3D * prmCL3Dprog
double getDoubleParam(const char *param, int arg=0)
double maxShiftX
Maximum shift X.
int Niter
Number of iterations.
double maxPsi
Maximum psi.
bool generateAlignedVolumes
Generate aligned subvolumes.
FileName fnSel
Input selfile with the images to quantify.
bool randomizeStartingOrientation
Use this option to avoid aligning at the beginning all the missing wedges.
double maxShiftY
Maximum shift Y.
int Ncodes
Final number of code vectors.
double DWTsparsity
DWT Sparsity factor (0<f<1; 1=drop all coefficients, 0=do not drop any coefficient) ...
const char * getParam(const char *param, int arg=0)
void max(Image< double > &op1, const Image< double > &op2)
int Nneighbours
Number of neighbours.
double maxShiftZ
Maximum shift Z.
bool dontAlign
Don&#39;t align.
void readParams(XmippProgram *program)
Definition: mask.cpp:1284
double maxTilt
Maximum tilt.
bool classifyAllImages
Clasify all images.
FileName fnOut
Output rootname.
bool checkParam(const char *param)
int Ncodes0
Initial number of code vectors.
double maxRot
Maximum rot.
FileName fnCodes0
Input selfile with initial codes.
int getIntParam(const char *param, int arg=0)
double PminSize
Minimum size of a node.
FileName fnSym
Symmetry file or code.
double maxFreq
Maximum frequency for the alignment.

◆ run()

void ProgClassifyCL3D::run ( )
virtual

Run.

Reimplemented from XmippProgram.

Definition at line 1483 of file mpi_classify_CLTomo_prog.cpp.

1484 {
1486  frmFunc = Python::getFunctionRef("sh_alignment.frm", "frm_align");
1487  wedgeClass = Python::getClassRef("sh_alignment.tompy.filter", "GeneralWedge");
1488 
1489  CREATE_LOG();
1490  show();
1491  produceSideInfo();
1492 
1493  // Run all iterations
1494  int level = 0;
1495  int Q = vq.P.size();
1496  vq.run(fnOut, level);
1497 
1498  while (Q < Ncodes)
1499  {
1500  if (node->rank == 0)
1501  std::cout << "Spliting nodes ...\n";
1502 
1503  int Nclean = vq.cleanEmptyNodes();
1504  int Nsplits = XMIPP_MIN(Q,Ncodes-Q) + Nclean;
1505 
1506  for (int i = 0; i < Nsplits; i++)
1507  {
1508  vq.splitFirstNode();
1509  if (node->rank == 0)
1510  std::cout << "Currently there are " << vq.P.size() << " nodes"
1511  << std::endl;
1512  }
1513 
1514  Q = vq.P.size();
1515  level++;
1516  vq.run(fnOut, level);
1517  }
1518  if (node->rank == 0)
1519  {
1520  std::sort(vq.P.begin(), vq.P.end(), SDescendingClusterSort());
1521  Q = vq.P.size();
1522  MetaDataDb SFq, SFclassified, SFaux, SFaux2;
1523  for (int q = 0; q < Q; q++)
1524  {
1525  SFq.read(formatString("class%06d_images@%s_classes_level_%02d.xmd", q + 1,
1526  fnOut.c_str(), level));
1527  SFq.fillConstant(MDL_REF, integerToString(q + 1));
1528  SFq.fillConstant(MDL_ENABLED, "1");
1529  SFclassified.unionAll(SFq);
1530  }
1531  SFaux = SF;
1532  SFaux.subtraction(SFclassified, MDL_IMAGE);
1533  SFaux.fillConstant(MDL_ENABLED, "-1");
1534  SFaux2.join1(SFclassified, SF, MDL_IMAGE, LEFT);
1535  SFclassified.clear();
1536  SFaux2.unionAll(SFaux);
1537  SFaux.clear();
1538  SFaux.sort(SFaux2, MDL_IMAGE);
1539  SFaux.write(fnOut + "_images.xmd");
1540  }
1541 
1542  // Produce aligned volumes
1544  {
1545  node->barrierWait();
1546  MetaDataDb MDimages, MDaligned;
1547  MDimages.read(fnOut+"_images.xmd");
1548  FileName fnAligned=fnOut+"_aligned.stk";
1549  FileName fnAlignedMD=fnOut+"_aligned.xmd";
1550  if (node->rank==0)
1551  createEmptyFile(fnAligned,(int)Xdim,(int)Ydim,(int)Zdim,MDimages.size());
1552  node->barrierWait();
1553  size_t idx=1;
1554  Image<double> V, Valigned;
1555  FileName fnImg;
1556  Matrix2D<double> A,E,T;
1557  Matrix1D<double> r(3);
1558  size_t itemId;
1559  for (size_t objId: MDimages.ids())
1560  {
1561  if (idx%node->size==node->rank)
1562  {
1563  MDimages.getValue(MDL_IMAGE,fnImg,objId);
1564  V.read(fnImg);
1565  V().setXmippOrigin();
1566  double x,y,z,rot,tilt,psi;
1567  MDimages.getValue(MDL_ANGLE_ROT,rot,objId);
1568  MDimages.getValue(MDL_ANGLE_TILT,tilt,objId);
1569  MDimages.getValue(MDL_ANGLE_PSI,psi,objId);
1570  MDimages.getValue(MDL_SHIFT_X,x,objId);
1571  MDimages.getValue(MDL_SHIFT_Y,y,objId);
1572  MDimages.getValue(MDL_SHIFT_Z,z,objId);
1573 
1574  Euler_angles2matrix(rot, tilt, psi, E, true);
1575  XX(r)=x;
1576  YY(r)=y;
1577  ZZ(r)=z;
1578  translation3DMatrix(r,T);
1579  A=E*T;
1580 
1581  applyGeometry(xmipp_transformation::BSPLINE3,Valigned(),V(),A,xmipp_transformation::IS_NOT_INV,xmipp_transformation::DONT_WRAP,0.);
1582  Valigned.write(fnAligned,idx,true,WRITE_REPLACE);
1583  }
1584  if (node->rank==0)
1585  {
1586  MDimages.getValue(MDL_ITEM_ID,itemId,objId);
1587  size_t newId=MDaligned.addObject();
1588  MDaligned.setValue(MDL_IMAGE,formatString("%06d@%s",idx,fnAligned.c_str()),newId);
1589  MDaligned.setValue(MDL_ITEM_ID,itemId,newId);
1590  }
1591  ++idx;
1592  }
1593 
1594  if (node->rank == 0)
1595  {
1596  MDaligned.write(fnAlignedMD);
1597 
1598  // Correct the class groups
1599  Q = vq.P.size();
1600  MetaDataDb SFq, SFqaligned;
1601  FileName fnClass;
1602  for (int q = 0; q < Q; q++)
1603  {
1604  fnClass=formatString("class%06d_images@%s_classes_level_%02d.xmd", q + 1, fnOut.c_str(), level);
1605  SFq.read(fnClass);
1607  SFqaligned.join1(SFq,MDaligned,MDL_ITEM_ID);
1608  SFqaligned.write(fnClass,MD_APPEND);
1609  }
1610  }
1611  }
1612  CLOSE_LOG();
1613 }
size_t size
Definition: xmipp_mpi.h:52
void subtraction(const MetaDataDb &mdIn, const MDLabel label)
Rotation angle of an image (double,degrees)
std::vector< size_t > objId
std::vector< CL3DClass * > P
List of nodes.
void Euler_angles2matrix(T alpha, T beta, T gamma, Matrix2D< T > &A, bool homogeneous)
Definition: geometry.cpp:624
PyObject * frmFunc
Pointer to the Python FRM alignment function.
bool getValue(MDObject &mdValueOut, size_t id) const override
Tilting angle of an image (double,degrees)
static double * y
#define CLOSE_LOG()
Shift for the image in the X axis (double)
void applyGeometry(int SplineDegree, MultidimArray< std::complex< double > > &V2, const MultidimArray< std::complex< double > > &V1, const Matrix2D< double > &A, bool inv, bool wrap, std::complex< double > outside, MultidimArray< double > *BcoeffsPtr)
#define CREATE_LOG()
void write(const FileName &name="", size_t select_img=ALL_IMAGES, bool isStack=false, int mode=WRITE_OVERWRITE, CastWriteMode castMode=CW_CAST, int _swapWrite=0)
Special label to be used when gathering MDs in MpiMetadataPrograms.
String integerToString(int I, int _width, char fill_with)
bool generateAlignedVolumes
Generate aligned subvolumes.
virtual IdIteratorProxy< false > ids()
void initPythonAndNumpy()
void barrierWait()
Definition: xmipp_mpi.cpp:171
doublereal * x
#define i
Is this image enabled? (int [-1 or 1])
Unique identifier for items inside a list or set (std::size_t)
void fillConstant(MDLabel label, const String &value) override
int Ncodes
Final number of code vectors.
Name of an image from which MDL_IMAGE is coming from.
#define XX(v)
Definition: matrix1d.h:85
template void translation3DMatrix(const Matrix1D< float > &translation, Matrix2D< float > &resMatrix, bool inverse)
size_t addObject() override
PyObject * getFunctionRef(const std::string &moduleName, const std::string &funcName)
void clear() override
Definition: metadata_db.cpp:54
double z
void createEmptyFile(const FileName &filename, int xdim, int ydim, int Zdim, size_t select_img, bool isStack, int mode, int _swapWrite, const MDRowVec *md)
void write(const FileName &outFile, WriteModeMetaData mode=MD_OVERWRITE) const override
void sort(struct DCEL_T *dcel)
Definition: sorting.cpp:18
size_t size() const override
#define XMIPP_MIN(x, y)
Definition: xmipp_macros.h:181
PyObject * getClassRef(const std::string &moduleName, const std::string &className)
size_t rank
Definition: xmipp_mpi.h:52
#define YY(v)
Definition: matrix1d.h:93
Class to which the image belongs (int)
void sort(MetaDataDb &MDin, const MDLabel sortLabel, bool asc=true, int limit=-1, int offset=0)
PyObject * wedgeClass
Pointer to the Python GeneralWedge class.
bool setValue(const MDObject &mdValueIn, size_t id) override
Definition: metadata_db.cpp:90
void run(const FileName &fnOut, int level)
double psi(const double x)
String formatString(const char *format,...)
Shift for the image in the Z axis (double)
void read(const FileName &inFile, const std::vector< MDLabel > *desiredLabels=NULL, bool decomposeStack=true) override
FileName fnOut
Output rootname.
void join1(const MetaDataDb &mdInLeft, const MetaDataDb &mdInRight, const MDLabel label, JoinType type=LEFT)
void unionAll(const MetaDataDb &mdIn)
int read(const FileName &name, DataMode datamode=DATA, size_t select_img=ALL_IMAGES, bool mapData=false, int mode=WRITE_READONLY)
Shift for the image in the Y axis (double)
void produceSideInfo()
Produce side info.
void renameColumn(MDLabel oldLabel, MDLabel newLabel) override
Name of an image (std::string)
#define ZZ(v)
Definition: matrix1d.h:101

◆ show()

void ProgClassifyCL3D::show ( ) const
virtual

Show.

Reimplemented from XmippProgram.

Definition at line 1347 of file mpi_classify_CLTomo_prog.cpp.

1348 {
1349  if (!verbose)
1350  return;
1351  std::cout << "Input images: " << fnSel << std::endl
1352  << "Output images: " << fnOut << std::endl
1353  << "Iterations: " << Niter << std::endl;
1354  if (fnCodes0!="")
1355  std::cout << "CodesSel0: " << fnCodes0 << std::endl;
1356  else
1357  std::cout << "Codes0: " << Ncodes0 << std::endl;
1358  std::cout << "Codes: " << Ncodes << std::endl
1359  << "Neighbours: " << Nneighbours << std::endl
1360  << "Minimum node size: " << PminSize << std::endl
1361  << "Sparsity: " << sparsity << std::endl
1362  << "DWT Sparsity: " << DWTsparsity << std::endl
1363  << "Maximum shift Z: " << maxShiftZ << std::endl
1364  << "Maximum shift Y: " << maxShiftY << std::endl
1365  << "Maximum shift X: " << maxShiftX << std::endl
1366  << "Maximum rot: " << maxRot << std::endl
1367  << "Maximum tilt: " << maxTilt << std::endl
1368  << "Maximum psi: " << maxPsi << std::endl
1369  << "Classify all images: " << classifyAllImages << std::endl
1370  << "Symmetry: " << fnSym << std::endl
1371  << "Don't align: " << dontAlign << std::endl
1372  << "Generate aligned volumes:" << generateAlignedVolumes << std::endl
1373  ;
1374  mask.show();
1375 }
double sparsity
Sparsity factor (0<f<1; 1=drop all coefficients, 0=do not drop any coefficient)
double maxShiftX
Maximum shift X.
int Niter
Number of iterations.
double maxPsi
Maximum psi.
bool generateAlignedVolumes
Generate aligned subvolumes.
FileName fnSel
Input selfile with the images to quantify.
double maxShiftY
Maximum shift Y.
int Ncodes
Final number of code vectors.
double DWTsparsity
DWT Sparsity factor (0<f<1; 1=drop all coefficients, 0=do not drop any coefficient) ...
int Nneighbours
Number of neighbours.
double maxShiftZ
Maximum shift Z.
bool dontAlign
Don&#39;t align.
int verbose
Verbosity level.
double maxTilt
Maximum tilt.
bool classifyAllImages
Clasify all images.
FileName fnOut
Output rootname.
int Ncodes0
Initial number of code vectors.
double maxRot
Maximum rot.
FileName fnCodes0
Input selfile with initial codes.
double PminSize
Minimum size of a node.
FileName fnSym
Symmetry file or code.
void show() const
Definition: mask.cpp:1021

Member Data Documentation

◆ classifyAllImages

bool ProgClassifyCL3D::classifyAllImages

Clasify all images.

Definition at line 259 of file mpi_classify_CLTomo.h.

◆ dontAlign

bool ProgClassifyCL3D::dontAlign

Don't align.

Definition at line 287 of file mpi_classify_CLTomo.h.

◆ DWTsparsity

double ProgClassifyCL3D::DWTsparsity

DWT Sparsity factor (0<f<1; 1=drop all coefficients, 0=do not drop any coefficient)

Definition at line 256 of file mpi_classify_CLTomo.h.

◆ fnCodes0

FileName ProgClassifyCL3D::fnCodes0

Input selfile with initial codes.

Definition at line 226 of file mpi_classify_CLTomo.h.

◆ fnOut

FileName ProgClassifyCL3D::fnOut

Output rootname.

Definition at line 229 of file mpi_classify_CLTomo.h.

◆ fnSel

FileName ProgClassifyCL3D::fnSel

Input selfile with the images to quantify.

Definition at line 223 of file mpi_classify_CLTomo.h.

◆ fnSym

FileName ProgClassifyCL3D::fnSym

Symmetry file or code.

Definition at line 232 of file mpi_classify_CLTomo.h.

◆ frmFunc

PyObject* ProgClassifyCL3D::frmFunc

Pointer to the Python FRM alignment function.

Definition at line 340 of file mpi_classify_CLTomo.h.

◆ gaussianInterpolator

GaussianInterpolator ProgClassifyCL3D::gaussianInterpolator

Definition at line 334 of file mpi_classify_CLTomo.h.

◆ generateAlignedVolumes

bool ProgClassifyCL3D::generateAlignedVolumes

Generate aligned subvolumes.

Definition at line 290 of file mpi_classify_CLTomo.h.

◆ mask

Mask ProgClassifyCL3D::mask

Mask.

Definition at line 284 of file mpi_classify_CLTomo.h.

◆ maxFreq

double ProgClassifyCL3D::maxFreq

Maximum frequency for the alignment.

Definition at line 250 of file mpi_classify_CLTomo.h.

◆ maxFreqMask

MultidimArray<unsigned char> ProgClassifyCL3D::maxFreqMask

MaxFreq mask.

Definition at line 349 of file mpi_classify_CLTomo.h.

◆ maxPsi

double ProgClassifyCL3D::maxPsi

Maximum psi.

Definition at line 281 of file mpi_classify_CLTomo.h.

◆ maxRot

double ProgClassifyCL3D::maxRot

Maximum rot.

Definition at line 275 of file mpi_classify_CLTomo.h.

◆ maxShift

double ProgClassifyCL3D::maxShift

Max shift.

Definition at line 346 of file mpi_classify_CLTomo.h.

◆ maxShiftX

double ProgClassifyCL3D::maxShiftX

Maximum shift X.

Definition at line 271 of file mpi_classify_CLTomo.h.

◆ maxShiftY

double ProgClassifyCL3D::maxShiftY

Maximum shift Y.

Definition at line 268 of file mpi_classify_CLTomo.h.

◆ maxShiftZ

double ProgClassifyCL3D::maxShiftZ

Maximum shift Z.

Definition at line 265 of file mpi_classify_CLTomo.h.

◆ maxTilt

double ProgClassifyCL3D::maxTilt

Maximum tilt.

Definition at line 278 of file mpi_classify_CLTomo.h.

◆ Ncodes

int ProgClassifyCL3D::Ncodes

Final number of code vectors.

Definition at line 241 of file mpi_classify_CLTomo.h.

◆ Ncodes0

int ProgClassifyCL3D::Ncodes0

Initial number of code vectors.

Definition at line 238 of file mpi_classify_CLTomo.h.

◆ Niter

int ProgClassifyCL3D::Niter

Number of iterations.

Definition at line 235 of file mpi_classify_CLTomo.h.

◆ Nneighbours

int ProgClassifyCL3D::Nneighbours

Number of neighbours.

Definition at line 244 of file mpi_classify_CLTomo.h.

◆ node

MpiNode* ProgClassifyCL3D::node

Definition at line 331 of file mpi_classify_CLTomo.h.

◆ objId

std::vector<size_t> ProgClassifyCL3D::objId

Definition at line 325 of file mpi_classify_CLTomo.h.

◆ PminSize

double ProgClassifyCL3D::PminSize

Minimum size of a node.

Definition at line 247 of file mpi_classify_CLTomo.h.

◆ randomizeStartingOrientation

bool ProgClassifyCL3D::randomizeStartingOrientation

Use this option to avoid aligning at the beginning all the missing wedges.

Definition at line 262 of file mpi_classify_CLTomo.h.

◆ SF

MetaDataDb ProgClassifyCL3D::SF

Definition at line 322 of file mpi_classify_CLTomo.h.

◆ SL

SymList ProgClassifyCL3D::SL

Definition at line 293 of file mpi_classify_CLTomo.h.

◆ sparsity

double ProgClassifyCL3D::sparsity

Sparsity factor (0<f<1; 1=drop all coefficients, 0=do not drop any coefficient)

Definition at line 253 of file mpi_classify_CLTomo.h.

◆ vq

CL3D ProgClassifyCL3D::vq

Definition at line 328 of file mpi_classify_CLTomo.h.

◆ wedgeClass

PyObject* ProgClassifyCL3D::wedgeClass

Pointer to the Python GeneralWedge class.

Definition at line 343 of file mpi_classify_CLTomo.h.

◆ Xdim

size_t ProgClassifyCL3D::Xdim

Definition at line 337 of file mpi_classify_CLTomo.h.

◆ Ydim

size_t ProgClassifyCL3D::Ydim

Definition at line 337 of file mpi_classify_CLTomo.h.

◆ Zdim

size_t ProgClassifyCL3D::Zdim

Definition at line 337 of file mpi_classify_CLTomo.h.


The documentation for this class was generated from the following files: