u1casc-ordinary/obs_phi2.hpp

   1 #ifndef OBS_PHI2_HPP
   2 #define OBS_PHI2_HPP
   3
   4 #include "latlib/o815/o815.h"
   5
   6 #include "latlib/writeout.h"
   7
   8 #include "latlib/obstat.hpp"
   9
  10 #include <iostream>
  11
  12 using namespace std;
  13
  14 class obs_phi2 : public o815::obs {
  15
  16 public:
  17   struct obsmem {
  18     double phi2[2];
  19   };
  20   obs_phi2(o815 *_O815);
  21
  22 private:
  23   void _start();
  24   void _meas(bool loadedobs, const int& nthmeas);
  25   void _finish();
  26
  27   obsmem* OM;
  28
  29   void phi2Compute();
  30   static double phi2Sus(vector< vector<double> > *vals, void *para);
  31
  32   sim *Sim;
  33   obstat<double,double> oPhi2[2];
  34 };
  35
  36 obs_phi2::obs_phi2(o815 *_O815) : o815::obs("phi2",
  37                                             _O815->paraQ->getParaNames() +
  38                                             "phi2one:phi2one_err:phi2onesus:phi2onesus_err"
  39                                             ":phi2two:phi2two_err:phi2twosus:phi2twosus_err",
  40                                             _O815, sizeof(obsmem) ) {
  41   OM = (obsmem*)obsMem;
  42   Sim = (sim*)O815->Sim;
  43 }
  44
  45 void obs_phi2::_start() {
  46   //*out << "OBS_test: start" << endl;
  47 };
  48
  49 void obs_phi2::_meas(bool loadedobs, const int& nthmeas) {
  50   if (!loadedobs)
  51     phi2Compute();
  52
  53   for (int iflav=0; iflav<2; iflav++)
  54     oPhi2[iflav].addMeas(OM->phi2[iflav]);
  55 };
  56
  57 void obs_phi2::_finish() {
  58   *out << O815->paraQ->getParaVals();
  59
  60   int compid_phi2[2], compid_phi2sus[2];
  61
  62   for (int iflav=0; iflav<2; iflav++) {
  63     compid_phi2[iflav] = oPhi2[iflav].computeEasy();
  64     compid_phi2sus[iflav] = oPhi2[iflav].computeJack(obs_phi2::phi2Sus, &(Sim->lsize4));
  65
  66     *out << "\t" << oPhi2[iflav].getMean(compid_phi2[iflav]) << "\t" << oPhi2[iflav].getErr(compid_phi2[iflav]);
  67     *out << "\t" << oPhi2[iflav].getMean(compid_phi2sus[iflav]) << "\t" << oPhi2[iflav].getErr(compid_phi2sus[iflav]);
  68
  69     oPhi2[iflav].reset();
  70   }
  71
  72   *out << endl;
  73 };
  74
  75 void obs_phi2::phi2Compute()
  76 {
  77   for (int iflav=0; iflav<2; iflav++) {
  78     OM->phi2[iflav] = 0;
  79
  80     for (int ix = 0; ix < Sim->lsize4; ix++)
  81       OM->phi2[iflav] += norm(Sim->phi[ iflav*Sim->lsize4 + ix ]);
  82
  83     OM->phi2[iflav] /= Sim->lsize4;
  84   }
  85 }
  86
  87 double obs_phi2::phi2Sus(vector< vector<double> > *vals, void *para) {
  88   double mean=0, mean2=0;
  89   for(vector< vector<double> >::iterator valIt = vals->begin(); valIt != vals->end(); ++valIt)
  90     {
  91       mean += (*valIt)[0];
  92       mean2 += pow((*valIt)[0],2);
  93     }
  94   mean /= vals->size();
  95   mean2 /= vals->size();
  96
  97   return ( mean2 - pow(mean,2) ) * *(int*)para;
  98 }
  99
 100 #endif