Something went wrong on our end
-
Bertram Kopf authored00f8a03a
JpsiGamKsKlKKHist.cc 11.62 KiB
#include <getopt.h>
#include <fstream>
#include <sstream>
#include <string>
#include "Examples/JpsiGamKsKlKK/JpsiGamKsKlKKHist.hh"
#include "Examples/JpsiGamKsKlKK/JpsiGamKsKlKKEventList.hh"
#include "Examples/JpsiGamKsKlKK/AbsJpsiGamKsKlKKLh.hh"
#include "TFile.h"
#include "TH1F.h"
#include "TH2F.h"
#include "TMath.h"
#include "TNtuple.h"
#include "ErrLogger/ErrLogger.hh"
JpsiGamKsKlKKHist::JpsiGamKsKlKKHist(boost::shared_ptr<const JpsiGamKsKlKKEventList> theEvtList) :
_theTFile(0),
_PhiPhiMassDataHist(0),
_PhiPhiMassMcHist(0),
_PhiPhiMassFittedHist(0),
_KpKmMassDataHist(0),
_KpKmMassMcHist(0),
_KpKmMassFittedHist(0),
_KsKlMassDataHist(0),
_KsKlMassMcHist(0),
_KsKlMassFittedHist(0),
_dataTuple(0),
_mcTuple(0)
{
if(0==theEvtList){
Alert <<"JpsiGamKsKlKKEventList* theEvtList is a 0 pointer !!!!" ; // << endmsg;
exit(1);
}
initRootStuff();
const std::vector<JpsiGamKsKlKKEvtData*> dataList=theEvtList->getDataVecs();
std::vector<JpsiGamKsKlKKEvtData*>::const_iterator it=dataList.begin();
while(it!=dataList.end())
{
plotDalitz(_dalitzDataHist, (*it), 1.);
plotPhiPhi(_PhiPhiMassDataHist, (*it), 1. );
plotKsKl( _KsKlMassDataHist, (*it), 1. );
plotKpKm( _KpKmMassDataHist, (*it), 1. );
plotCostPhiKs( _costKs_KsKlHeliDataHist, _phiKs_KsKlHeliDataHist,(*it), 1. );
plotCostPhiKs( _costKp_KpKmHeliDataHist, _phiKp_KpKmHeliDataHist,(*it), 1. );
plotCostGam( _costGamCmDataHist,(*it), 1. );
fillTuple(_dataTuple, (*it), 1.);
++it;
}
const std::vector<JpsiGamKsKlKKEvtData*> mcList=theEvtList->getMcVecs();
it=mcList.begin();
while(it!=mcList.end())
{
plotDalitz(_dalitzMcHist, (*it), 1.);
plotDalitz(_dalitzMcHist, (*it), 1.);
plotPhiPhi(_PhiPhiMassMcHist, (*it), 1. );
plotKsKl( _KsKlMassMcHist, (*it), 1. );
plotKpKm( _KpKmMassMcHist, (*it), 1. );
plotCostPhiKs( _costKs_KsKlHeliMcHist, _phiKs_KsKlHeliMcHist,(*it), 1. );
plotCostPhiKp( _costKp_KpKmHeliMcHist, _phiKp_KpKmHeliMcHist, (*it), 1. );
plotCostGam( _costGamCmMcHist, (*it), 1. );
fillTuple(_mcTuple, (*it), 1.);
++it;
}}
JpsiGamKsKlKKHist::JpsiGamKsKlKKHist(boost::shared_ptr<AbsJpsiGamKsKlKKLh> theJpsiGamKsKlKKLh, paramGamKsKlKK& fitParam) :
_theTFile(0),
_dalitzDataHist(0),
_dalitzMcHist(0),
_dalitzFittedHist(0),
_PhiPhiMassDataHist(0),
_PhiPhiMassMcHist(0),
_PhiPhiMassFittedHist(0),
_KpKmMassDataHist(0),
_KpKmMassMcHist(0),
_KpKmMassFittedHist(0),
_KsKlMassDataHist(0),
_KsKlMassMcHist(0),
_KsKlMassFittedHist(0),
_dataTuple(0),
_mcTuple(0)
{
if(0==theJpsiGamKsKlKKLh){
Alert <<"JpsiGamKsKlKKLh* theJpsiGamKsKlKKLh is a 0 pointer !!!!" ; // << endmsg;
exit(1);
}
initRootStuff();
boost::shared_ptr<const JpsiGamKsKlKKEventList> theEvtList=theJpsiGamKsKlKKLh->getEventList();
const std::vector<JpsiGamKsKlKKEvtData*> dataList=theEvtList->getDataVecs();
std::vector<JpsiGamKsKlKKEvtData*>::const_iterator it=dataList.begin();
while(it!=dataList.end())
{
plotDalitz(_dalitzDataHist, (*it), 1.);
plotPhiPhi(_PhiPhiMassDataHist, (*it), 1. );
plotKsKl( _KsKlMassDataHist, (*it), 1. );
plotKpKm( _KpKmMassDataHist, (*it), 1. );
plotCostPhiKs( _costKs_KsKlHeliDataHist, _phiKs_KsKlHeliDataHist, (*it), 1. );
plotCostPhiKp( _costKp_KpKmHeliDataHist, _phiKp_KpKmHeliDataHist, (*it), 1. );
plotCostGam( _costGamCmDataHist, (*it), 1. );
fillTuple(_dataTuple, (*it), 1.);
++it;
}
const std::vector<JpsiGamKsKlKKEvtData*> mcList=theEvtList->getMcVecs();
it=mcList.begin();
while(it!=mcList.end())
{
plotDalitz(_dalitzMcHist, (*it), 1.);
plotPhiPhi(_PhiPhiMassMcHist, (*it), 1. );
plotKsKl( _KsKlMassMcHist, (*it), 1. );
plotKpKm( _KpKmMassMcHist, (*it), 1. );
plotCostPhiKs( _costKs_KsKlHeliMcHist, _phiKs_KsKlHeliMcHist, (*it), 1. );
plotCostPhiKp( _costKp_KpKmHeliMcHist, _phiKp_KpKmHeliMcHist,(*it), 1. );
plotCostGam( _costGamCmMcHist,(*it), 1. );
double evtWeight= theJpsiGamKsKlKKLh->calcEvtIntensity((*it), fitParam);
plotDalitz(_dalitzFittedHist, (*it),evtWeight );
plotPhiPhi(_PhiPhiMassFittedHist, (*it), evtWeight );
plotKsKl( _KsKlMassFittedHist, (*it), evtWeight );
plotKpKm( _KpKmMassFittedHist, (*it), evtWeight );
plotCostPhiKs( _costKs_KsKlHeliFittedHist, _phiKs_KsKlHeliFittedHist,(*it), evtWeight );
plotCostPhiKp( _costKp_KpKmHeliFittedHist, _phiKp_KpKmHeliFittedHist,(*it), evtWeight );
plotCostGam( _costGamCmFittedHist,(*it), evtWeight );
fillTuple(_mcTuple, (*it), evtWeight);
++it;
}
double integralData=(double) theEvtList->getDataVecs().size();
Info <<"No of fit data events " << integralData ; // << endmsg;
double integralFitted=(double) theEvtList->getMcVecs().size();
Info <<"No of fit events " << integralFitted ; // << endmsg;
Info <<"scaling factor " << integralData/integralFitted ; // << endmsg;
double scaleFactor = integralData/integralFitted;
_dalitzFittedHist->Scale(scaleFactor);
_PhiPhiMassFittedHist->Scale(scaleFactor);
_KsKlMassFittedHist->Scale(scaleFactor);
_KpKmMassFittedHist->Scale(scaleFactor);
_costKs_KsKlHeliFittedHist->Scale(scaleFactor);
_phiKs_KsKlHeliFittedHist->Scale(scaleFactor);
_costKp_KpKmHeliFittedHist->Scale(scaleFactor);
_phiKp_KpKmHeliFittedHist->Scale(scaleFactor);
_costGamCmFittedHist->Scale(scaleFactor);
}
JpsiGamKsKlKKHist::~JpsiGamKsKlKKHist()
{
_theTFile->Write();
_theTFile->Close();
}
void JpsiGamKsKlKKHist::initRootStuff()
{
std::string rootFileName="./JpsiGamKsKlKK.root";
_theTFile=new TFile(rootFileName.c_str(),"recreate");
int xbins = 50;
double xmin=0.8;
double xmax =10.;
int ybins=xbins;
double ymin=xmin;
double ymax=xmax;
_dalitzDataHist= new TH2F("_dalitzDataHist","Dpl K+K- K+#pi^{0} data",xbins, xmin, xmax, ybins, ymin, ymax );
_dalitzMcHist= new TH2F("_dalitzMcHist","Dpl K+K- K+#pi^{0} MC",xbins, xmin, xmax, ybins, ymin, ymax);
_dalitzFittedHist= new TH2F("_dalitzFittedHist","Dpl K+K- K+#pi^{0} fit",xbins, xmin, xmax, ybins, ymin, ymax );
int nbins=50;
xmin=2;
xmax=3.2;
_PhiPhiMassDataHist = new TH1F("_PhiPhiMassDataHist", "_PhiPhiMassDataHist", nbins, xmin, xmax);
_PhiPhiMassMcHist = new TH1F("_PhiPhiMassMcHist", "_PhiPhiMassMcHist", nbins, xmin, xmax);
_PhiPhiMassFittedHist = new TH1F("_PhiPhiMassFittedHist", "_PhiPhiMassFittedHist", nbins, xmin, xmax);
xmin=0.8;
xmax=1.2;
_KsKlMassDataHist = new TH1F("_KsKlMassDataHist", "_KsKlMassDataHist", nbins, xmin, xmax);
_KsKlMassMcHist = new TH1F("_KsKlMassMcHist", "_KsKlMassMcHist", nbins, xmin, xmax);
_KsKlMassFittedHist = new TH1F("_KsKlMassFittedHist", "_KsKlMassFittedHist", nbins, xmin, xmax);
_KpKmMassDataHist = new TH1F("_KpKmMassDataHist", "_KpKmMassDataHist", nbins, xmin, xmax);
_KpKmMassMcHist = new TH1F("_KpKmMassMcHist", "_KpKmMassMcHist", nbins, xmin, xmax);
_KpKmMassFittedHist = new TH1F("_KpKmMassFittedHist", "_KpKmMassFittedHist", nbins, xmin, xmax);
_costKs_KsKlHeliDataHist= new TH1F("_costKs_KsKlHeliDataHist", "cos(#Theta_{Ks}) KsKlHeli data", 100, -1., 1.);
_costKs_KsKlHeliMcHist= new TH1F("_costKs_KsKlHeliMcHist", "cos(#Theta_{Ks}) KsKlHeli Mc", 100, -1., 1.);
_costKs_KsKlHeliFittedHist= new TH1F("_costKs_KsKlHeliFittedHist", "cos(#Theta_{Ks}) KsKlHeli Fit", 100, -1, 1);
_phiKs_KsKlHeliDataHist= new TH1F("_phiKs_KsKlHeliDataHist", "#Phi_{Ks} KsKlHeli data", 100, -1., 1.);
_phiKs_KsKlHeliMcHist= new TH1F("_phiKs_KsKlHeliMcHist", "#Phi_{Ks} KsKlHeli Mc", 100, -1., 1.);
_phiKs_KsKlHeliFittedHist= new TH1F("_phiKs_KsKlHeliFittedHist", "#Phi_{Ks} KsKlHeli Fit", 100, -1, 1);
_costKp_KpKmHeliDataHist= new TH1F("_costKp_KpKmHeliDataHist", "cos(#Theta_{K+}) KsKlHeli data", 100, -1., 1.);
_costKp_KpKmHeliMcHist= new TH1F("_costKp_KpKmHeliMcHist", "cos(#Theta_{K+}) KsKlHeli Mc", 100, -1., 1.);
_costKp_KpKmHeliFittedHist= new TH1F("_costKp_KpKmHeliFittedHist", "cos(#Theta_{K+}) KsKlHeli Fit", 100, -1, 1);
_phiKp_KpKmHeliDataHist= new TH1F("_phiKp_KpKmHeliDataHist", "#Phi_{K+} KsKlHeli data", 100, -1., 1.);
_phiKp_KpKmHeliMcHist= new TH1F("_phiKp_KpKmHeliMcHist", "#Phi_{K+} KsKlHeli Mc", 100, -1., 1.);
_phiKp_KpKmHeliFittedHist= new TH1F("_phiKp_KpKmHeliFittedHist", "#Phi_{K+} KsKlHeli Fit", 100, -1, 1);
_costGamCmDataHist= new TH1F("_costGamCmDataHist", "cos(#Theta_{#gamma}) CM data", 100, -1., 1.);
_costGamCmMcHist= new TH1F("_costGamCmMcHist", "cos(#Theta_{#gamma}) CM Mc", 100, -1., 1.);
_costGamCmFittedHist= new TH1F("_costGamCmFittedHist", "cos(#Theta_{#gamma}) CM Fit", 100, -1, 1);
std::string tupleVariables = "mKsKlKpKm:mKsKl:mKpKm:KsKlKpKmCostTheta:gamCosTheta:KsKlCosTheta:KpKmCosTheta:KsCosTheta:KpCosTheta:weight";
_dataTuple=new TNtuple("_dataTuple", "data ntuple", tupleVariables.data());
_mcTuple=new TNtuple("_mcTuple", "mc ntuple", tupleVariables.data());
}
void JpsiGamKsKlKKHist::plotDalitz(TH2F* theHisto, const JpsiGamKsKlKKEvtData* theData, double weight)
{
Vector4<float> V4_KsKlKpKm_HeliPsi;
Vector4<float> V4_KsKl_HeliPsi = theData->V4_KsKl_HeliPsi ;
Vector4<float> V4_KpKm_HeliPsi = theData->V4_KpKm_HeliPsi ;
Vector4<float> V4_gamma_HeliPsi = theData->V4_gamma_HeliPsi ;
double gphi1 = (V4_gamma_HeliPsi + V4_KsKl_HeliPsi ).M();
double gphi2 = (V4_gamma_HeliPsi + V4_KpKm_HeliPsi ).M();
theHisto->Fill( gphi1*gphi1, gphi2*gphi2 ,weight);
}
void JpsiGamKsKlKKHist::plotPhiPhi(TH1F* theHisto, const JpsiGamKsKlKKEvtData* theData, double weight){
Vector4<float> v4 = theData->V4_KsKlKpKm_HeliPsi;
theHisto->Fill( v4.M(), weight );
}
void JpsiGamKsKlKKHist::plotKsKl(TH1F* theHisto, const JpsiGamKsKlKKEvtData* theData, double weight){
Vector4<float> v4 = theData->V4_KsKl_HeliPsi;
theHisto->Fill( v4.M(), weight );
}
void JpsiGamKsKlKKHist::plotKpKm(TH1F* theHisto, const JpsiGamKsKlKKEvtData* theData, double weight){
Vector4<float> v4 = theData->V4_Kp_HeliKpKm;
theHisto->Fill( v4.M(), weight );
}
void JpsiGamKsKlKKHist::plotCostPhiKs(TH1F* theCostHisto, TH1F* thePhiHisto, const JpsiGamKsKlKKEvtData* theData, double weight){
Vector4<float> v4 = theData->V4_Ks_HeliKsKl;
theCostHisto->Fill( v4.CosTheta(), weight );
thePhiHisto->Fill( v4.Phi(), weight );
}
void JpsiGamKsKlKKHist::plotCostPhiKp(TH1F* theCostHisto, TH1F* thePhiHisto, const JpsiGamKsKlKKEvtData* theData, double weight){
Vector4<float> v4 = theData->V4_Kp_HeliKpKm;
theCostHisto->Fill( v4.CosTheta(), weight );
thePhiHisto->Fill( v4.Phi(), weight );
}
void JpsiGamKsKlKKHist::plotCostGam(TH1F* theCostHisto, const JpsiGamKsKlKKEvtData* theData, double weight){
Vector4<float> v4 = theData->V4_gamma_HeliPsi;
theCostHisto->Fill( v4.CosTheta(), weight );
}
void JpsiGamKsKlKKHist::fillTuple( TNtuple* theTuple, const JpsiGamKsKlKKEvtData* theData, double weight)
{
Vector4<float> V4_KsKlKpKm_HeliPsi = theData->V4_KsKlKpKm_HeliPsi ;
Vector4<float> V4_KsKl_HeliPsi= theData->V4_KsKl_HeliPsi ;
Vector4<float> V4_KpKm_HeliPsi= theData->V4_KpKm_HeliPsi ;
Vector4<float> V4_gamma_HeliPsi= theData->V4_gamma_HeliPsi ;
Vector4<float> V4_KsKl_HeliKsKlKpKm= theData->V4_KsKl_HeliKsKlKpKm ;
Vector4<float> V4_KpKm_HeliKsKlKpKm= theData->V4_KpKm_HeliKsKlKpKm ;
Vector4<float> V4_Ks_HeliKsKl= theData->V4_Ks_HeliKsKl ;
Vector4<float> V4_Kp_HeliKpKm= theData->V4_Kp_HeliKpKm ;
theTuple->Fill( V4_KsKlKpKm_HeliPsi.M(),
V4_KsKl_HeliPsi.M(),
V4_KpKm_HeliPsi.M(),
V4_KsKlKpKm_HeliPsi.CosTheta(),
V4_gamma_HeliPsi.CosTheta(),
V4_KsKl_HeliKsKlKpKm.CosTheta(),
V4_KpKm_HeliKsKlKpKm.CosTheta(),
V4_Ks_HeliKsKl.CosTheta(),
V4_Kp_HeliKpKm.CosTheta(),
weight);
}