PolnProfileFitAnalysis.h

//-*-C++-*-
/***************************************************************************
 *
 *   Copyright (C) 2005 by Willem van Straten
 *   Licensed under the Academic Free License version 2.1
 *
 ***************************************************************************/
 
// psrchive/More/Polarimetry/Pulsar/PolnProfileFitAnalysis.h
 
#ifndef __Pulsar_PolnProfileFitAnalysis_h
#define __Pulsar_PolnProfileFitAnalysis_h
 
#include "Pulsar/PolnProfileFit.h"
#include "MEAL/StokesError.h"
#include "MEAL/StokesCovariance.h"
 
#include "MEAL/MuellerTransformation.h"
 
#include <vector>
 
namespace Pulsar {

  class PolnProfileFit::Analysis {
 
  public:

    static bool verbose;

    Analysis ();

    void set_compute_error (bool flag = true);

    void set_fit (PolnProfileFit*);

    Estimate<double> get_relative_error () const;

    Estimate<double> get_multiple_correlation () const;

    Estimate<double> get_relative_conditional_error () const;

    double get_c_varphi ();

    double get_C_varphi ();

    double get_expected_relative_error (std::vector<unsigned>& histogram);

    double get_c_varphi_error () const;

    double get_Rmult ();

    double get_cond_var ();

    Matrix<8,8,double> get_correlation () const;

    void output_C_varphi (const char* filename);
 
  protected:

    Estimate<double> relative_error;

    Estimate<double> relative_conditional_error;

    Estimate<double> multiple_correlation;

    bool built;

    void build (bool only_relative_error = false);

    Reference::To<PolnProfileFit> fit;

    Reference::To<MEAL::Complex2> xform;

    unsigned max_harmonic;

    Jones<double> del_deleta (unsigned i, const Jones<double>& K) const;

    double delR2_varphiJ_delS (Matrix<8,8,double>& delC_delS);

    void delC_delS( Matrix<8,8,double>& delC_delSre,
                    Matrix<8,8,double>& delC_delSim, unsigned k ) const;

    Jones<double> delrho_delS (unsigned k) const;

    void get_curvature (Matrix<8,8,double>& curvature);
    void add_curvature (Matrix<8,8,double>& curvature);

    void set_harmonic (unsigned index);

    void evaluate (unsigned index);
 
    Jones<double> model_result;
    std::vector< Jones<double> > model_gradient;
 
    Jones<double> xform_result;
    std::vector< Jones<double> > xform_gradient;
 
    Jones<double> phase_result;
    std::vector< Jones<double> > phase_gradient;

    MEAL::StokesError error;

    MEAL::StokesCovariance stokes_covariance;

    Matrix<8,8,double> covariance;

    Matrix <7,7,double> C_JJ;

    Matrix <7,7,double> inv_C_JJ;

    Vector <7,double> C_varphiJ;
    
    double c_varphi;

    double R2_varphiJ;
 
    double efac;

    bool compute_error;

    std::vector< Reference::To<MEAL::Complex2> > inputs;
    
    Matrix<4,4,double> Xi;
 
  private:
 
    void initialize();
 
    std::vector< double > weights;
    std::vector< double > store_covariance;
 
    void compute_weights (unsigned nharmonic);
 
  };
 

  class ScalarProfileFitAnalysis {
 
  public:

    void set_fit (const PolnProfileFit*);

    void set_spectrum (const Profile*, double rescale_variance = 1.0);

    void set_max_harmonic (unsigned max_harmonic);

    void set_variance (double v);

    double get_mean_variance () const
    { return 0.5 * (variance.real() + variance.imag()); }

    std::complex<double> get_variance () const { return variance; }

    Estimate<double> get_error () const;

    void get_curvature (Matrix<2,2,double>& curvature);

    void delC_delS ( Matrix<2,2,double>& delC_delSre,
                     Matrix<2,2,double>& delC_delSim,
                     unsigned index ) const;
 
  protected:

    Reference::To<const Profile> spectrum;

    const std::complex<float>* amps;

    unsigned max_harmonic;
 
    Matrix<2,2,double> covariance;
    std::complex<double> variance;
  };
 
}
 
#endif