ml_align2d.h

Go to the documentation of this file.

/***************************************************************************
 *
 * Authors:    Sjors Scheres           scheres@cnb.csic.es (2004)
 *
 * Unidad de  Bioinformatica of Centro Nacional de Biotecnologia , CSIC
 *
 * This program is free software; you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation; either version 2 of the License, or
 * (at your option) any later version.
 *
 * This program is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License
 * along with this program; if not, write to the Free Software
 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA
 * 02111-1307  USA
 *
 *  All comments concerning this program package may be sent to the
 *  e-mail address 'xmipp@cnb.csic.es'
 ***************************************************************************/

#ifndef _MLALIGN2D_H
#define _MLALIGN2D_H

#include "ml2d.h"

#define FOR_ALL_THREAD_REFNO() \
int refno, load; \
while ((load = getThreadRefnoJob(refno)) > 0) \
    for (int i = 0; i < load; i++, refno = (refno + 1) % model.n_ref)
#define FOR_ALL_THREAD_REFNO_NODECL() \
while ((load = getThreadRefnoJob(refno)) > 0) \
    for (int i = 0; i < load; i++, refno = (refno + 1) % model.n_ref)

class ProgML2D;
// This structure is needed to pass parameters to the threads
typedef struct
{
    int thread_id;
    ProgML2D * prm;
}
structThreadTasks;

#define SPECIAL_ITER 0

void * doThreadsTasks(void * data);

class ProgML2D: public ML2DBaseProgram
{

public:
  bool no_iem;

    MultidimArray<int> mask, omask;
    int threadTask;
    barrier_t barrier, barrier2, barrier3;
    pthread_t * th_ids;
    structThreadTasks * threads_d;

    double LL, sumfracweight;
    double wsum_sigma_noise, wsum_sigma_offset, sumw_allrefs;
    std::vector<double> sumw, sumw2, sumwsc, sumwsc2, sumw_mirror;
    std::vector<MultidimArray<double > > wsum_Mref;
    MultidimArray<int> Msignificant;
    MultidimArray<double> Mimg;
    std::vector<double> allref_offsets;
    std::vector<double> pdf_directions;
    std::vector<MultidimArray<double> > mref;
    std::vector<MultidimArray<std::complex<double> > > fref;
    std::vector<MultidimArray<std::complex<double > > > wsumimgs;
    int opt_refno, iopt_psi, iopt_flip;
    double trymindiff, opt_scale, bgmean, opt_psi;
    double fracweight, maxweight2, dLL;


    std::vector<MultidimArray<std::complex<double> > > Fimg_flip, mysumimgs;
    std::vector<double> refw, refw2, refwsc2, refw_mirror, sumw_refpsi;
    double wsum_corr, sum_refw, sum_refw2, maxweight;
    double wsum_sc, wsum_sc2, wsum_offset, old_bgmean;
    double mindiff;
    size_t sigdim;
    int ioptx, iopty;
    std::vector<int> ioptx_ref, iopty_ref, ioptflip_ref;
    double pfs_mindiff;
    MultidimArray<double> pfs_maxweight;
    MultidimArray<double> pfs_weight;
    int pfs_count;

    std::vector<double> A2;
    double Xi2;

    //These are for refno work assigns to threads
    size_t refno_index, refno_load_param;
    int refno_load, refno_count;
    // Which group does this image belong to in iteration 0 (generation of K references)
    int mygroup;
    void readParams();
    void defineParams();

public:
    ProgML2D();
    virtual void show();
    virtual void produceSideInfo();
    virtual void produceSideInfo2();

    void calculatePdfInplane();

    void rotateReference();

    void reverseRotateReference();

    void preselectLimitedDirections(double &phi, double &theta);

    void preselectFastSignificant();

    void expectationSingleImage(Matrix1D<double> &opt_offsets);

    /*** Threads functions */
    void createThreads();

    void destroyThreads();

    int getThreadRefnoJob(int &refno);

    void awakeThreads(ThreadTask task, int start_refno, int load = 1);

    void doThreadRotateReferenceRefno();

    void doThreadReverseRotateReferenceRefno();

    void doThreadPreselectFastSignificantRefno();

    void doThreadExpectationSingleImageRefno();

    void doThreadESIUpdateRefno();

    virtual void iteration();

    virtual void expectation();

    void maximizeModel(ModelML2D &model);

    virtual void maximization();

    void correctScaleAverage();

    virtual void addPartialDocfileData(const MultidimArray<double> &data, size_t first, size_t last);

    virtual void writeOutputFiles(const ModelML2D &model, OutputType outputType = OUT_FINAL);

    virtual void readModel(ModelML2D &model, int block);
    FileName getBaseName(String suffix = "", int number = -1);

    virtual void printModel(const String &msg, const ModelML2D & model);

};



#endif