u1casc-ordinary/obs_corrphichimass.hpp

   1 #ifndef OBS_CORRPHICHIMASS_HPP
   2 #define OBS_CORRPHICHIMASS_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 #include <complex>
  12
  13 #include <cmath>
  14
  15 #define U1 0
  16 #define U2 1
  17 #define F1 2
  18 #define F2 3
  19
  20 using namespace std;
  21
  22 class obs_corrphichimass : public o815::obs {
  23
  24 public:
  25   obs_corrphichimass(o815 *_O815);
  26
  27 private:
  28   void _start();
  29   void _meas(bool loadedobs, const int& nthmeas);
  30   void _finish();
  31
  32   void corrCompute();
  33   static complex<double> effMass(vector < complex<double> > *preCalculated, vector< complex<double> > *excludedMeas, int nmeas, void *para);
  34   static void preEffMass(vector< vector < complex<double> > > *allVals, vector < complex<double> > *preCalculated, void *para);
  35
  36   sim *Sim;
  37   obstat< complex<double>, complex<double>  > oC;
  38
  39   int spatialV;
  40
  41   complex<double> *OM;
  42 };
  43
  44 obs_corrphichimass::obs_corrphichimass(o815 *_O815) : o815::obs("corrphichimass",
  45                                                                   _O815->paraQ->getParaNames() +
  46                                                                   "tsep"
  47                                                                   ":phichi_effmass:phichi_effmass_err",
  48                                                                   _O815, sizeof(complex<double>)*(_O815->comargs.lsize[1]/2+1) ) {
  49
  50   OM  = (complex<double>*)(obsMem);
  51
  52   Sim = (sim*)O815->Sim;
  53   spatialV = O815->comargs.lsize[0] * O815->comargs.lsize[0] * O815->comargs.lsize[0];
  54 }
  55
  56 void obs_corrphichimass::_start() {
  57   //*out << "OBS_test: start" << endl;
  58 };
  59
  60 void obs_corrphichimass::_meas(bool loadedobs, const int& nthmeas) {
  61   if (!loadedobs)
  62     corrCompute();
  63
  64   oC.addMeas( OM, O815->comargs.lsize[1]/2+1 );
  65 };
  66
  67 void obs_corrphichimass::_finish() {
  68   int compid_corr[O815->comargs.lsize[1]/2][4];
  69   int compid_effmass[O815->comargs.lsize[1]/2-1][4];
  70
  71   for (int itsep = 0; itsep < O815->comargs.lsize[1]/2; itsep++)
  72     compid_corr[itsep][U1] = oC.computeEasy(itsep);
  73
  74   for (int itsep = 0; itsep < O815->comargs.lsize[1]/2-1; itsep++) {
  75     pair<int,int> effmasspass( itsep, O815->comargs.lsize[1]/2 );
  76     compid_effmass[itsep][U1] = oC.computeJack(obs_corrphichimass::preEffMass, obs_corrphichimass::effMass, &effmasspass);
  77   }
  78
  79   for (int itsep = 0; itsep < O815->comargs.lsize[1]/2-1; itsep++) {
  80     *out << O815->paraQ->getParaVals();
  81     *out << "\t" << itsep;
  82
  83     *out << "\t" << real( oC.getMean(compid_effmass[itsep][U1]) )
  84          << "\t" << real( oC.getErr (compid_effmass[itsep][U1]) );
  85
  86     *out << endl;
  87   }
  88
  89   oC.reset();
  90 };
  91
  92 void obs_corrphichimass::corrCompute()
  93 {
  94   complex<double> phislice[O815->comargs.lsize[1]];
  95
  96   OM[O815->comargs.lsize[1]/2] = 0;
  97
  98   for (int it = 0; it < O815->comargs.lsize[1]; it++) {
  99     phislice[it] = 0;
 100
 101     for (int ix = 0; ix < spatialV; ix++)
 102       phislice[it] += Sim->phi[ 0*Sim->lsize4 + it*spatialV + ix ] * Sim->phi[ 1*Sim->lsize4 + it*spatialV + ix ];
 103       //phislice[it] += conj(Sim->phi[ 0*Sim->lsize4 + it*spatialV + ix ]) * Sim->phi[ 1*Sim->lsize4 + it*spatialV + ix ];
 104
 105     phislice[it] /= spatialV;
 106
 107     OM[O815->comargs.lsize[1]/2] += phislice[it];
 108   }
 109
 110   for (int itsep = 0; itsep < O815->comargs.lsize[1]/2; itsep++) {
 111     OM[itsep] = 0;
 112
 113     for (int it = 0; it < O815->comargs.lsize[1]; it++)
 114       OM[itsep] += phislice[ (it+itsep)%O815->comargs.lsize[1] ] * conj( phislice[it] );
 115
 116     OM[itsep] /= O815->comargs.lsize[1];
 117   }
 118
 119   OM[O815->comargs.lsize[1]/2] /= O815->comargs.lsize[1];
 120 }
 121
 122 void obs_corrphichimass::preEffMass(vector< vector < complex<double> > > *allVals, vector < complex<double> > *preCalculated, void *para) {
 123   pair<int,int> *myparas = (pair<int,int>*)para;
 124
 125   preCalculated->push_back(0);
 126   preCalculated->push_back(0);
 127   preCalculated->push_back(0);
 128
 129   for(vector< vector< complex<double> > >::iterator valIt = allVals->begin(); valIt != allVals->end(); ++valIt) {
 130     (*preCalculated)[0] += (*valIt)[myparas->first];
 131     (*preCalculated)[1] += (*valIt)[myparas->first+1];
 132     (*preCalculated)[2] += (*valIt)[myparas->second];
 133   }
 134 }
 135
 136 complex<double> obs_corrphichimass::effMass(vector < complex<double> > *preCalculated, vector< complex<double> > *excludedMeas, int nmeas, void *para) {
 137   pair<int,int> *myparas = (pair<int,int>*)para;
 138
 139   double disconnected = norm( ( (*preCalculated)[2] - (*excludedMeas)[myparas->second] ) / (complex<double>)(nmeas-1) );
 140
 141   return std::log( abs(
 142                        ( ( (*preCalculated)[0] - (*excludedMeas)[myparas->first] ) / (complex<double>)(nmeas-1) - disconnected ) /
 143                        ( ( (*preCalculated)[1] - (*excludedMeas)[myparas->first+1] ) / (complex<double>)(nmeas-1) - disconnected )
 144                         ) );
 145 }
 146
 147 #endif