ProfileShiftFit.h

//-*-C++-*-
/***************************************************************************
 *
 *   Copyright (C) 2008 by Paul Demorest
 *   Licensed under the Academic Free License version 2.1
 *
 ***************************************************************************/
 
#ifndef __Pulsar_ProfileShiftFit_h
#define __Pulsar_ProfileShiftFit_h
 
#include "Pulsar/Algorithm.h"
#include "FTransform.h"
#include "Estimate.h"
#include "toa.h"
#include "BoxMuller.h"
 
namespace Pulsar
{
  class Profile;
  class Integration;
 
 
  class ProfileShiftFit : public Algorithm
  {   
 
  public:
 
    enum Error_Method {
      Traditional_Chi2,     // Chi-squared 2nd derivatives
      MCMC_Variance,        // Shift variance from MCMC
      Numerical,            // Shift variance from numerical integration
    };
 
    ProfileShiftFit ();
 
    ~ProfileShiftFit();
 
    void reset();
 
    void set_nharm(unsigned nh);
 
    unsigned get_nharm();
 
    void set_error_method(Error_Method e) { err_meth = e; }
 
    void set_mcmc_iterations(int nit) { mcmc_it = nit; }
 
    int get_mcmc_iterations() const { return mcmc_it; }
 
    void set_standard (const Profile* p);
 
    void set_Profile (const Profile* p);
 
    void compute();
 
    Tempo::toa toa (const Integration*);
 
    Estimate<double> get_shift();
 
    Estimate<double> get_scale();
 
    double get_mse();
 
    double get_sigma2();
 
    double get_snr();
 
    double get_reduced_chisq () const;
 
 
    double get_effective_duty_cycle () const;
 
    void apply_scale_and_shift(Profile *p);
 
    int mcmc_trials;
    int mcmc_accept;
 
    bool choose_maximum_harmonic;
 
  protected:
 
    void init();
 
    void choose_nharm ();
 
    unsigned nharm;
 
    unsigned effective_nharm;
 
    Reference::To<const Profile> std;
 
    unsigned nbins_std;
 
    float *fstd;
 
    double std_pow;
 
    Reference::To<const Profile> prof;
 
    unsigned nbins_prof;
 
    float *fprof;
 
    float *fccf;
 
    unsigned nbins_ccf;
 
    
    double ccf(double phi);
 
    double dccf(double phi);
 
    double d2ccf(double phi);
 
    double log_shift_pdf(double phi);
 
    double log_shift_pdf_pos(double phi);
 
    Error_Method err_meth;
 
    void error_traditional();
 
    void error_numerical();
 
    void error_mcmc_pdf_var();
 
    int mcmc_it;
 
    bool computed;
 
    double shift;
 
    double eshift;
 
    double correction;
 
    double scale;
 
    double escale;
 
    double sigma2;
 
    double mse;
 
    double chi2;
 
    unsigned dof;
 
    double snr;
 
  private:
 
    BoxMuller nrand;
 
    double mcmc_state;
 
    double mcmc_log_pdf;
 
    void mcmc_init();
 
    double mcmc_sample();
 
    double max_log_pdf;
 
    static double f_pdf(double phi, void *_psf);
    static double f_pdf_x2(double phi, void *_psf);
 
  };
 
}
 
 
#endif