Xmipp  v3.23.11-Nereus
Public Member Functions | Public Attributes | List of all members
ProgAlignTiltPairs Class Reference
Libraries » Reconstruction library » align_tilt_pairs (Align tilted and untilted images in a random conical tilt experiment)

#include <align_tilt_pairs.h>

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

Public Member Functions

void defineParams ()
 Define parameters in the command line. More...
 
void readParams ()
 Read parameters from the command line. More...
 
void show ()
 Show. More...
 
void processImage ()
 Process a single image. More...
 
void run ()
 Run over the whole input metadata. More...
 
bool centerTiltedImage (const MultidimArray< double > &imgU, bool flip, double inPlaneU, double shiftXu, double shiftYu, double alphaT, double alphaU, double tilt, const MultidimArray< double > &imgT, double &shiftX, double &shiftY, CorrelationAux &auxCorr)
 Center one tilted image. 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)
 
ParamDef * getParamDef (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 show () 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 fnIn
 
FileName fnRef
 
FileName fnOut
 
double max_shift
 
bool do_stretch
 
bool do_not_align_tilted
 
- 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
 
ProgramDef * progDef
 Program definition and arguments parser. More...
 
std::map< String, CommentList > defaultComments
 
int argc
 Original command line arguments. More...
 
const char ** argv
 

Detailed Description

Centilt parameters.

Definition at line 41 of file align_tilt_pairs.h.

Member Function Documentation

◆ centerTiltedImage()

bool ProgAlignTiltPairs::centerTiltedImage ( const MultidimArray< double > &  imgU,
bool  flip,
double  inPlaneU,
double  shiftXu,
double  shiftYu,
double  alphaT,
double  alphaU,
double  tilt,
const MultidimArray< double > &  imgT,
double &  shiftX,
double &  shiftY,
CorrelationAux &  auxCorr 
)

Center one tilted image.

Definition at line 90 of file align_tilt_pairs.cpp.

96 {
97  // Cosine stretching, store stretched image in imgT2DClass
98  MultidimArray<double> imgT2DClass;
99  Matrix2D<double> E, E2D, Mu2D, A2D;
100 
101  if (!do_stretch)
102  tilt = flip ? 180. : 0.;
103 
104  Euler_angles2matrix(flip ? alphaU : -alphaU, tilt, alphaT, E, true);
105 
106  rotation2DMatrix(inPlaneU, Mu2D, true);
107  MAT_ELEM(Mu2D,0,2) = flip ? shiftXu : -shiftXu;
108  MAT_ELEM(Mu2D,1,2) = -shiftYu;
109 
110  if (flip)
111  {
112  // We need to multiply the first column by -1, because
113  // in Xmipp the Mu2D order is: first flip (if necessary), then rot, then shift
114  //MAT_ELEM(Mu2D,0,0)*=-1;
115 
116  MAT_ELEM(Mu2D,1,0)*=-1;
117  MAT_ELEM(Mu2D,2,0)*=-1;
118 
119  // We need to multiply the first row by -1, as needed by the RCT theory
120  MAT_ELEM(Mu2D,0,1)*=-1;
121  MAT_ELEM(Mu2D,0,2)*=-1;
122  }
123 
124  E2D.initIdentity(3);
125  MAT_ELEM(E2D,0,0)=MAT_ELEM(E,0,0);
126  MAT_ELEM(E2D,0,1)=MAT_ELEM(E,0,1);
127  MAT_ELEM(E2D,1,0)=MAT_ELEM(E,1,0);
128  MAT_ELEM(E2D,1,1)=MAT_ELEM(E,1,1);
129  A2D = Mu2D * E2D.inv();
130 
131  applyGeometry(xmipp_transformation::LINEAR, imgT2DClass, imgT, A2D, xmipp_transformation::IS_NOT_INV, xmipp_transformation::WRAP);
132 
133  // Calculate best shift
134  auto max_shift_pixels=(int)(max_shift/100.0*XSIZE(imgT));
135  CorrelationAux aux;
136  double corr=bestShift(imgRefU, imgT2DClass, shiftX, shiftY, auxCorr, nullptr, max_shift_pixels);
137  // double corr=bestShiftRealSpace(imgRefU, imgT2DClass, shiftX, shiftY, NULL, max_shift_pixels, 0.5);
138 
139 #ifdef DEBUG
140 
141  std::cout << "alphaU=" << alphaU << " inplaneU=" << inPlaneU << " tilt=" << tilt << " alphaT=" << alphaT << std::endl;
142  std::cout << "Best shift= " << shiftX << " " << shiftY << " corr=" << corr << std::endl;
143 #endif
144 
145  Matrix1D<double> vShift(2);
146  XX(vShift)=shiftX;
147  YY(vShift)=shiftY;
148  Matrix2D<double> Tt, Tt2D;
149  translation2DMatrix(vShift, Tt2D,true);
150  Tt=A2D.inv()*Tt2D.inv()*A2D;
151 
152  // Readjust shift
153  shiftX=MAT_ELEM(Tt,0,2);
154  shiftY=MAT_ELEM(Tt,1,2);
155  double shift = sqrt(shiftX * shiftX + shiftY * shiftY);
156 
157 //#define DEBUG_2
158 #ifdef DEBUG_2
159 
160  Image<double> save;
161  save() = imgT;
162  save.write("PPPtilted.xmp");
163  save() = imgRefU;
164  save.write("PPPuntiltedRef.xmp");
165  save() = imgT2DClass;
166  save.write("PPPtiltedAdjusted.xmp");
167  std::cout << "Corrected shift= " << shiftX << " " << shiftY << std::endl;
168  std::cout << "Press any key\n";
169  char c = getchar();
170 #endif
171 
172  return (shift < (double)max_shift_pixels) || corr<0;
173 }
translation2DMatrix
template void translation2DMatrix(const Matrix1D< float > &, Matrix2D< float > &, bool inverse)
Euler_angles2matrix
void Euler_angles2matrix(T alpha, T beta, T gamma, Matrix2D< T > &A, bool homogeneous)
Definition: geometry.cpp:624
c
doublereal * c
Definition: numerical_recipes.cpp:2230
sqrt
void sqrt(Image< double > &op)
Definition: image_operate.cpp:210
xmipp_transformation::LINEAR
Definition: transformations_defines.h:32
applyGeometry
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)
Definition: transformations.cpp:668
ImageBase::write
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)
Definition: xmipp_image_base.cpp:292
Matrix2D::inv
void inv(Matrix2D< T > &result) const
Definition: matrix2d.cpp:663
MAT_ELEM
#define MAT_ELEM(m, i, j)
Definition: matrix2d.h:116
XX
#define XX(v)
Definition: matrix1d.h:85
XSIZE
#define XSIZE(v)
Definition: multidim_array_base.h:91
rotation2DMatrix
void rotation2DMatrix(T ang, Matrix2D< T > &result, bool homogeneous)
Definition: transformations.cpp:361
YY
#define YY(v)
Definition: matrix1d.h:93
bestShift
template float bestShift(MultidimArray< float > &, float &, float &, const MultidimArray< int > *, int)
Matrix2D::initIdentity
void initIdentity()
Definition: matrix2d.h:673

◆ defineParams()

void ProgAlignTiltPairs::defineParams ( )
virtual

Define parameters in the command line.

Reimplemented from XmippProgram.

Definition at line 37 of file align_tilt_pairs.cpp.

38 {
39  //Usage
40  addUsageLine("Center the tilted images of all tilted-untilted image pairs.");
41  addUsageLine("+This program receives as input a metadata with two sets of images, untilted and tilted.");
42  addUsageLine("+Untilted images must have been previously classified and aligned. Then, the tilted images");
43  addUsageLine("+are aligned to their untilted versions. The alignment parameters from the classification plus");
44  addUsageLine("+the alignment parameters from the tilted-untilted comparison are used to deduce the 3D");
45  addUsageLine("+alignment parameters for the tilted images.");
46 
47  //Examples
48  addExampleLine("To center tilted images allowing a maximum shift of 10 pixels:",false);
49  addExampleLine("xmipp_align_tilt_pairs -i class_00001@classes.xmd -o alignedImages.xmd --max_shift 10");
50 
51  // Params
52  addParamsLine(" -i <metadata> : Input metadata with untilted and tilted images");
53  addParamsLine(" -o <metadata> : Output metadata file with rotations & translations for 3D reconstruction.");
54  addParamsLine(" --ref <file> : 2D average of the untilted images");
55  addParamsLine(" [--max_shift <value=10>] : Discard images shifting more than a given threshold (in percentage of the image size).");
56  addParamsLine(" :+Set it to 0 for no shift estimate between tilted and untilted images");
57  addParamsLine(" [--do_stretch] : Stretch tilted image to fit into the untilted one.");
58  addParamsLine(" :+Do it only with thin particles");
59  addParamsLine(" [--do_not_align_tilted] : Do not align tilted images to untilted ones.");
60  addParamsLine(" :+Do not align if the quality of the tilted images is low.");
61 }
XmippProgram::addExampleLine
void addExampleLine(const char *example, bool verbatim=true)
Definition: xmipp_program.cpp:286
XmippProgram::addUsageLine
void addUsageLine(const char *line, bool verbatim=false)
Definition: xmipp_program.cpp:282
XmippProgram::addParamsLine
void addParamsLine(const String &line)
Definition: xmipp_program.cpp:305

◆ processImage()

void ProgAlignTiltPairs::processImage ( )

Process a single image.

◆ readParams()

void ProgAlignTiltPairs::readParams ( )
virtual

Read parameters from the command line.

Reimplemented from XmippProgram.

Definition at line 64 of file align_tilt_pairs.cpp.

65 {
66  fnIn = getParam("-i");
67  fnOut = getParam("-o");
68  fnRef = getParam("--ref");
69  max_shift = getDoubleParam("--max_shift");
70  do_stretch = checkParam("--do_stretch");
71  do_not_align_tilted = checkParam("--do_not_align_tilted");
72 }
XmippProgram::getDoubleParam
double getDoubleParam(const char *param, int arg=0)
Definition: xmipp_program.cpp:341
XmippProgram::getParam
const char * getParam(const char *param, int arg=0)
Definition: xmipp_program.cpp:321
XmippProgram::checkParam
bool checkParam(const char *param)
Definition: xmipp_program.cpp:379

◆ run()

void ProgAlignTiltPairs::run ( )
virtual

Run over the whole input metadata.

Reimplemented from XmippProgram.

Definition at line 177 of file align_tilt_pairs.cpp.

178 {
179  Image<double> imgT;
180  MultidimArray<double> Maux;
181  Matrix2D<double> A(3, 3);
182 
183  MetaDataVec mdIn, mdOut;
184  mdIn.read(fnIn);
185  mdIn.removeDisabled();
186  if (!mdIn.containsLabel(MDL_ANGLE_TILT))
187  REPORT_ERROR(ERR_ARG_INCORRECT,"Input metadata does not have tilt information");
188 
189  initProgress(mdIn.size());
190 
191  Image<double> imgRef;
192  imgRef.read(fnRef);
193  imgRef().setXmippOrigin();
194 
195  size_t imgno = 0;
196  FileName fnTilted, fnUntilted;
197 
198  double alphaU, alphaT, tilt, inPlaneU, minusInPlaneU;
199  size_t nDiscarded = 0;
200  Matrix2D<double> E, Tu, Tup;
201  Matrix1D<double> vShift(3);
202  vShift.initZeros();
203  CorrelationAux auxCorr;
204  bool flip;
205 
206  for (const auto& row : mdIn)
207  {
208  row.getValue(MDL_FLIP, flip);
209  // Read untilted and tilted images
210  row.getValue(MDL_ANGLE_PSI, inPlaneU);
211  row.getValue(MDL_IMAGE, fnUntilted);
212  row.getValue(MDL_IMAGE_TILTED, fnTilted);
213  imgT.read(fnTilted);
214  imgT().setXmippOrigin();
215 
216 
217 #ifdef DEBUG
218  {
219  Image<double> imgU;
220  imgU.read(fnUntilted);
221  imgU().setXmippOrigin();
222  Image<double> save;
223  save() = imgU();
224  save.write("PPPuntilted.xmp");
225  Matrix2D<double> E;
226  save.readApplyGeo(fnUntilted, row);
227  save.write("PPPuntiltedAligned.xmp");
228  geo2TransformationMatrix(row, E);
229  //rotation2DMatrix(-inPlaneU,E);
230  //save() = imgU();
231  //applyGeometry(LINEAR, save(), imgU(), E, IS_NOT_INV, WRAP);
232  //save.write("PPPuntiltedAligned2.xmp");
233  }
234 #endif
235 
236  // Get alignment parameters
237  double shiftXu, shiftYu;
238  row.getValue(MDL_ANGLE_Y, alphaU);
239  row.getValue(MDL_ANGLE_Y2, alphaT);
240  row.getValue(MDL_ANGLE_TILT, tilt);
241  row.getValue(MDL_SHIFT_X, shiftXu);
242  row.getValue(MDL_SHIFT_Y, shiftYu);
243 
244  if (flip)
245  {
246  tilt += 180;
247  minusInPlaneU = inPlaneU+alphaU;
248  }
249  else
250  minusInPlaneU=-(inPlaneU+alphaU);
251 
252  // Correct untilted alignment
253  Euler_angles2matrix( minusInPlaneU, tilt, alphaT, E, true);
254  XX(vShift) = shiftXu;
255  YY(vShift) = shiftYu;
256  translation3DMatrix(vShift, Tu);
257  Tup = E * Tu * E.inv();
258 
259  // Align tilt and untilted projections
260  double shiftX = 0, shiftY = 0;
261  bool enable = true;
262  if (max_shift > 0)
263  enable = centerTiltedImage(imgRef(), flip, inPlaneU, shiftXu, shiftYu,
264  alphaT, alphaU, tilt, imgT(), shiftX, shiftY, auxCorr);
265  if (!enable)
266  {
267  ++nDiscarded;
268  shiftX = shiftY = 0;
269  }
270 
271  // Write results
272  MDRowVec rowOut;
273  rowOut.setValue(MDL_IMAGE, fnTilted);
274  rowOut.setValue(MDL_ANGLE_ROT, minusInPlaneU);
275  rowOut.setValue(MDL_ANGLE_TILT, tilt);
276  rowOut.setValue(MDL_ANGLE_PSI, alphaT);
277  rowOut.setValue(MDL_SHIFT_X, -MAT_ELEM(Tup,0,3) + shiftX);
278  rowOut.setValue(MDL_SHIFT_Y, -MAT_ELEM(Tup,1,3) + shiftY);
279  rowOut.setValue(MDL_ENABLED, (int) enable);
280  mdOut.addRow(rowOut);
281 
282  setProgress(++imgno);
283  }
284 
285  if (verbose && nDiscarded)
286  std::cout << " Discarded " << nDiscarded
287  << " tilted images that shifted too much" << std::endl;
288 
289  // Write out selfile
290  mdOut.write(fnOut);
291 }
MDL_ANGLE_ROT
Rotation angle of an image (double,degrees)
Definition: metadata_label.h:55
ProgAlignTiltPairs::centerTiltedImage
bool centerTiltedImage(const MultidimArray< double > &imgU, bool flip, double inPlaneU, double shiftXu, double shiftYu, double alphaT, double alphaU, double tilt, const MultidimArray< double > &imgT, double &shiftX, double &shiftY, CorrelationAux &auxCorr)
Center one tilted image.
Definition: align_tilt_pairs.cpp:90
MetaDataVec::read
void read(const FileName &inFile, const std::vector< MDLabel > *desiredLabels=nullptr, bool decomposeStack=true) override
Definition: metadata_vec.cpp:120
REPORT_ERROR
#define REPORT_ERROR(nerr, ErrormMsg)
Definition: xmipp_error.h:211
geo2TransformationMatrix
void geo2TransformationMatrix(const MDRow &imageGeo, Matrix2D< double > &A, bool only_apply_shifts)
Definition: transformations.cpp:178
Euler_angles2matrix
void Euler_angles2matrix(T alpha, T beta, T gamma, Matrix2D< T > &A, bool homogeneous)
Definition: geometry.cpp:624
MDL_ANGLE_TILT
Tilting angle of an image (double,degrees)
Definition: metadata_label.h:59
MDRowVec::setValue
void setValue(const MDObject &object) override
Definition: metadata_row_vec.cpp:199
MDL_SHIFT_X
Shift for the image in the X axis (double)
Definition: metadata_label.h:431
XmippProgram::initProgress
void initProgress(size_t total, size_t stepBin=60)
Definition: xmipp_program.cpp:251
ImageBase::write
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)
Definition: xmipp_image_base.cpp:292
Matrix2D::inv
void inv(Matrix2D< T > &result) const
Definition: matrix2d.cpp:663
MDL_ANGLE_PSI
Special label to be used when gathering MDs in MpiMetadataPrograms.
Definition: metadata_label.h:51
MetaDataVec::write
void write(const FileName &outFile, WriteModeMetaData mode=MD_OVERWRITE) const
Definition: metadata_vec.cpp:140
MDL_IMAGE_TILTED
Name of the tilted images associated to MDL_IMAGE.
Definition: metadata_label.h:252
MDL_ANGLE_Y2
Angle between y-axis and tilt-axis (double, degrees) for tilted micrographs.
Definition: metadata_label.h:67
ImageBase::readApplyGeo
int readApplyGeo(const FileName &name, const MDRow &row, const ApplyGeoParams &params=DefaultApplyGeoParams)
Definition: xmipp_image_base.cpp:264
MetaDataVec::size
size_t size() const override
Definition: metadata_vec.cpp:393
MDL_ENABLED
Is this image enabled? (int [-1 or 1])
Definition: metadata_label.h:231
MetaDataVec::addRow
size_t addRow(const MDRow &row) override
Definition: metadata_vec.cpp:225
MAT_ELEM
#define MAT_ELEM(m, i, j)
Definition: matrix2d.h:116
XX
#define XX(v)
Definition: matrix1d.h:85
translation3DMatrix
template void translation3DMatrix(const Matrix1D< float > &translation, Matrix2D< float > &resMatrix, bool inverse)
ERR_ARG_INCORRECT
Incorrect argument received.
Definition: xmipp_error.h:113
MDL_ANGLE_Y
Angle between y-axis and tilt-axis (double, degrees) for untilted micrographs.
Definition: metadata_label.h:66
MDL_FLIP
Flip the image? (bool)
Definition: metadata_label.h:236
XmippProgram::verbose
int verbose
Verbosity level.
Definition: xmipp_program.h:161
YY
#define YY(v)
Definition: matrix1d.h:93
XmippProgram::setProgress
void setProgress(size_t value=0)
Definition: xmipp_program.cpp:263
MetaData::removeDisabled
virtual void removeDisabled()
Definition: metadata_base.cpp:123
ImageBase::read
int read(const FileName &name, DataMode datamode=DATA, size_t select_img=ALL_IMAGES, bool mapData=false, int mode=WRITE_READONLY)
Definition: xmipp_image_base.cpp:119
MDL_SHIFT_Y
Shift for the image in the Y axis (double)
Definition: metadata_label.h:435
MetaDataVec::containsLabel
bool containsLabel(const MDLabel label) const override
Definition: metadata_vec.cpp:397
MDL_IMAGE
Name of an image (std::string)
Definition: metadata_label.h:246

◆ show()

void ProgAlignTiltPairs::show ( )

Show.

Definition at line 75 of file align_tilt_pairs.cpp.

76 {
77  std::cout
78  << "Input metadata: " << fnIn << std::endl
79  << "Output metadata: " << fnOut << std::endl
80  << "Reference image: " << fnRef << std::endl;
81  if (max_shift != 0)
82  std::cout
83  << "Discard images that shift more than: " << max_shift << std::endl;
84  if (do_stretch)
85  std::cout << "Stretching tilted images\n";
86 }

Member Data Documentation

◆ do_not_align_tilted

bool ProgAlignTiltPairs::do_not_align_tilted

Do not align tilted

Definition at line 55 of file align_tilt_pairs.h.

◆ do_stretch

bool ProgAlignTiltPairs::do_stretch

Do stretch

Definition at line 53 of file align_tilt_pairs.h.

◆ fnIn

FileName ProgAlignTiltPairs::fnIn

Filename input document file

Definition at line 45 of file align_tilt_pairs.h.

◆ fnOut

FileName ProgAlignTiltPairs::fnOut

Filename output document file

Definition at line 49 of file align_tilt_pairs.h.

◆ fnRef

FileName ProgAlignTiltPairs::fnRef

Filename untilted average

Definition at line 47 of file align_tilt_pairs.h.

◆ max_shift

double ProgAlignTiltPairs::max_shift

Discard images that shift more than max_shift

Definition at line 51 of file align_tilt_pairs.h.

The documentation for this class was generated from the following files:
Generated by   doxygen 1.8.13