// $Id$ // // File: JEventProcessor_strangeness.cc // Created: Tue Aug 29 08:59:06 EDT 2017 // Creator: staylor (on Linux ifarm1402.jlab.org 3.10.0-327.el7.x86_64 x86_64) // #include "JEventProcessor_strangeness.h" using namespace jana; #include #include #include "../include/mass_cuts.h" // Routine used to create our JEventProcessor #include #include extern "C"{ void InitPlugin(JApplication *app){ InitJANAPlugin(app); app->AddProcessor(new JEventProcessor_strangeness()); } } // "C" #include //------------------ // JEventProcessor_strangeness (Constructor) //------------------ JEventProcessor_strangeness::JEventProcessor_strangeness() { } //------------------ // ~JEventProcessor_strangeness (Destructor) //------------------ JEventProcessor_strangeness::~JEventProcessor_strangeness() { } //------------------ // init //------------------ jerror_t JEventProcessor_strangeness::init(void) { japp->CreateLock("mykklock"); BETA_CUT=0.9; // This is called once at program startup. EMULATE_TRIGGER=true; gPARMS->SetDefaultParameter("STRANGENESS:EMULATE_TRIGGER",EMULATE_TRIGGER); DROP_ONE_CELL_BCAL_SHOWERS=true; gPARMS->SetDefaultParameter("STRANGENESS:DROP_ONE_CELL_BCAL_SHOWERS",DROP_ONE_CELL_BCAL_SHOWERS); VETO_GAP_EVENTS=true; gPARMS->SetDefaultParameter("STRANGENESS:VETO_GAP_EVENTS",VETO_GAP_EVENTS); PROTON_GAMMA_DT_CUT=2.; gPARMS->SetDefaultParameter("STRANGENESS:PROTON_GAMMA_DT_CUT", PROTON_GAMMA_DT_CUT); SPLIT_CUT=0.5; gPARMS->SetDefaultParameter("STRANGENESS:SPLIT_CUT",SPLIT_CUT); PI0_VETO_CUT=0.042; // 3 sigma? gPARMS->SetDefaultParameter("STRANGENESS:PI0_VETO_CUT",PI0_VETO_CUT); FCAL_RADIAL_CUT=104.0; gPARMS->SetDefaultParameter("STRANGENESS:FCAL_RADIAL_CUT",FCAL_RADIAL_CUT); BCAL_Z_CUT=384.0; gPARMS->SetDefaultParameter("STRANGENESS:BCAL_Z_CUT",BCAL_Z_CUT); FILL_ROO_DATASET=false; gPARMS->SetDefaultParameter("STRANGENESS:FILL_ROO_DATASET", FILL_ROO_DATASET); NUM_SIGMA_BG=2; gPARMS->SetDefaultParameter("STRANGENESS:NUM_SIGMA_BG",NUM_SIGMA_BG); PI0_CUT_VALUE=0.015; ETA_CUT_VALUE=0.02; ETAPRIME_CUT_VALUE=0.02; gPARMS->SetDefaultParameter("STRANGENESS:PI0_CUT_VALUE",PI0_CUT_VALUE); gPARMS->SetDefaultParameter("STRANGENESS:ETA_CUT_VALUE",ETA_CUT_VALUE); gPARMS->SetDefaultParameter("STRANGENESS:ETAPRIME_CUT_VALUE",ETAPRIME_CUT_VALUE); gDirectory->mkdir("strangeness")->cd(); /* if (FILL_ROO_DATASET){ RooMass=new RooRealVar("m","",0.,3.,"GeV"); RooT=new RooRealVar("t","",0.,3.,"GeV^{-2}"); RooE=new RooRealVar("E","",3.,12.,"GeV"); RooWeight=new RooRealVar("weight","",-1.,1.,""); RooData=new RooDataSet("MassData","MassData",RooArgSet(*RooMass,*RooT,*RooE,*RooWeight)); } */ TrigBits = new TH1F("TrigBits","Trigger bits",10,-0.5,9.5); gDirectory->mkdir("MC")->cd(); { MCProtonThetaVsP=new TH2F("MCProtonThetaVsP","Thrown #theta vs p for proton",200,0,2,90,0,90); MCProtonThetaVsP->Sumw2(); MCProtonThetaVsP->SetXTitle("p [GeV/c]"); MCProtonThetaVsP->SetYTitle("#theta [degrees]"); MCProtonAcceptedP=new TH1F("MCProtonAcceptedP","Accepted proton",200,0,2); MCProtonAcceptedP->Sumw2(); MCBeam = new TH1F("MCBeam","Photon energy",48,2.8,12.4); MCBeamTagged = new TH1F("MCBeamTagged","Photon energy",48,2.8,12.4); // MCProtonAcceptance=new TH1F("MCProtonAcceptance","proton acceptance", //200,0,2); MCProtonP=new TH1F("MCProtonP","thrown proton distribution", 200,0,2); MCProtonP->Sumw2(); MCMissingMass =new TH1F("MCMissingMass","mX",1600,0,4.); MCMissingMassVsEbeam=new TH2F("MCMissingMassVsEbeam", "MX vs E#gamma", 48,2.8,12.4,1600,0,4.); for (unsigned int i=0;i<48;i++){ char myname[80]; sprintf(myname,"KpCosThetaH_vs_t_%3.1fGeV_mc",2.9+0.2*i); KpCosThetaH_vs_t_mc[i]=new TH2F(myname,myname,50,0,2,100,-1,1); sprintf(myname,"KpCosThetaVsPhi_mc_t%d",i); KpCosThetaVsPhi_mc[i]=new TH2F(myname,"cos(#theta) vs #phi",50,-M_PI,M_PI,50,-1,1); } } gDirectory->cd("../"); gDirectory->mkdir("PID")->cd(); { ProtonGammaTimeDiff=new TH1F("ProtonGammaTimeDiff","#deltat between proton and #gamma at vertex",200,-10,10); UnusedEnergy=new TH1F("UnusedEnergy","unused shower energy",100,0,0.5); NeutralBeta=new TH1F("NeutralBeta","#beta for neutrals",100,0,1.5); KinFitCL=new TH1F("KinFitCL","fit probability",100,0,1); KinFitChiSq=new TH1F("KinFitChiSq","#chi^2/ndf",1000,0,100); } gDirectory->cd("../"); gDirectory->mkdir("KpKmP")->cd(); { KpKmVertex=new TH2F("KpKmVertex","vertex",1000,0,100,200,0,20); KpKmReducedChiSq=new TH1F("KpKmReducedChiSq","#chi^2/ndf",1000,0,100); KpThetaVsPhi_lab=new TH2F("KpThetaVsPhi_lab","#theta vs #phi", 360,-180,180,180,0,90); KpKmMass= new TH1F("KpKmMass","K+K- mass",800,0.9,2.9); KpKmMass->SetXTitle("M(K^{+}K^{-}) [GeV]"); KpKmMass->SetYTitle("counts / 0.0025 GeV"); KmPMass=new TH1F("KmPMass","pK- mass",1600,0,4.); KmPMass->SetXTitle("M(pK^{-}) [GeV]"); KmPMass->SetYTitle("counts / 0.0025 GeV"); KpKmPDalitz=new TH2F("KpKmPDalitz","pK+K- Dalitz plot",400,0.9,8.9,700,2,18); KpKmMass_vs_Ebeam=new TH2F("KpKmMass_vs_Ebeam","M(K+K-) vs E#gamma", 48,2.8,12.4,500,0.95,1.45); PKm_vs_t=new TH2F("PKm_vs_t","M(pK-) vs t",50,0,2,200,1.4,2.4); gDirectory->mkdir("t-dependence")->cd(); for (unsigned int i=0;i<48;i++){ char myname[80]; sprintf(myname,"KpKm_vs_t_%3.1fGeV_v2",2.9+0.2*i); KpKm_vs_t[i]=new TH2F(myname,myname,50,0,2,300,0.9,1.2); sprintf(myname,"KpCosThetaH_vs_t_%3.1fGeV",2.9+0.2*i); KpCosThetaH_vs_t[i]=new TH2F(myname,myname,50,0,2,100,-1,1); sprintf(myname,"KpCosThetaVsPhi_t%d",i); KpCosThetaVsPhi[i]=new TH2F(myname,"cos(#theta) vs #phi",50,-M_PI,M_PI,50,-1,1); } gDirectory->cd("../"); } gDirectory->cd("../"); gDirectory->mkdir("KpKm1Gamma")->cd(); { KpKmGammaMass= new TH1F("KpKmGammaMass","K+K-#gamma mass",1600,0,4.); PhiGammaMass= new TH1F("PhiGammaMass","#phi#gamma mass",1600,0,4.); KpKmMass_1g= new TH1F("KpKmMass_1g","K+K- mass",800,0.9,2.9); KpKmMass_1g->SetXTitle("M(K^{+}K^{-}) [GeV]"); KmPMass_1g=new TH1F("KmPMass_1g","pK- mass",1600,0,4.); KpKmGammaMass_vs_Ebeam=new TH2F("KpKmGammaMass_vs_Ebeam", "M(K+K-#gamma) vs E#gamma", 48,2.8,12.4,1600,0,4.); } gDirectory->cd("../"); gDirectory->mkdir("KpKmP2Gamma")->cd(); { KpKm2GammaMass= new TH1F("KpKm2GammaMass","K+K-2#gamma mass",1600,0,4.); TwoGammaMass_KpKm= new TH1F("TwoGammaMass_KpKm","2#gamma mass",1600,0,4.); KpKmMass_2g= new TH1F("KpKmMass_2g","K+K- mass",800,0.9,2.9); KpKmMass_2g->SetXTitle("M(K^{+}K^{-}) [GeV]"); KpKmMass_2g->SetYTitle("counts / 0.0025 GeV"); PhiPi0Mass= new TH1F("PhiPi0Mass","K+K-2#gamma mass",1600,0,4.); PhiPi0Mass->SetXTitle("M(K^{+}K^{-})#pi^{0} [GeV]"); PhiPi0Mass->SetYTitle("counts / 0.0025 GeV"); PhiEtaMass= new TH1F("PhiEtaMass","K+K-2#gamma mass",1600,0,4.); PhiEtaMass->SetXTitle("M(K^{+}K^{-}#eta) [GeV]"); PhiEtaMass->SetYTitle("counts / 0.0025 GeV"); PhiEtaPrimeMass= new TH1F("PhiEtaPrimeMass","K+K-2#gamma mass",1600,0,4.); PhiEtaPrimeMass->SetXTitle("M(K^{+}K^{-}#eta') [GeV]"); PhiEtaPrimeMass->SetYTitle("counts / 0.0025 GeV"); KpKmPi0Mass= new TH1F("KpKmPi0Mass","K+K-2#gamma mass",1600,0,4.); KpKmPi0Mass->SetXTitle("M(K^{+}K^{-}#pi^{0}) [GeV]"); KpKmPi0Mass->SetYTitle("counts / 0.0025 GeV"); KpKmPi0Mass_baryon_cut= new TH1F("KpKmPi0Mass_baryon_cut","K+K-2#gamma mass",1600,0,4.); KpKmPi0Mass_baryon_cut->SetXTitle("M(K^{+}K^{-}#pi^{0}) [GeV]"); KpKmPi0Mass_baryon_cut->SetYTitle("counts / 0.0025 GeV"); KstarKm= new TH1F("KstarKm","K+K-2#gamma mass",1600,0,4.); KstarKp= new TH1F("KstarKp","K+K-2#gamma mass",1600,0,4.); KpKmEtaMass= new TH1F("KpKmEtaMass","K+K-2#gamma mass",1600,0,4.); KpKmEtaMass->SetXTitle("M(K^{+}K^{-}#eta) [GeV]"); KpKmEtaMass->SetYTitle("counts / 0.0025 GeV"); KpKmPi0Dalitz=new TH2F("KpKmPi0Dalitz","M^{2}(K+#pi0) vs M^{2}(K-#pi0)", 200,0,4,200,0,4); KpKm_with_Pi0= new TH1F("KpKm_with_Pi0","K+K- mass",800,0.9,2.9); KpKmMass_vs_2GammaMass=new TH2F("KpKmMass_vs_2GammaMass","K+K- mass vs 2#gamma mass",300,0.,1.5,600,0.9,2.5); PPi0_vs_KpKm=new TH2F("PPi0_vs_KpKm","p#pi^{0} mass vs K+K- mass",800,0.,4.,800,0.,4.); } gDirectory->cd("../"); gDirectory->mkdir("KpKmP3Gamma")->cd(); { KpKm3GammaMass= new TH1F("KpKm3GammaMass","K+K-3#gamma mass",1600,0,4.); TwoGammaMass_kpkm3g =new TH1F("TwoGammaMass_kpkm3g","M(2#gamma)",1600,0,4.); ThreeGammaMass_KpKm= new TH1F("ThreeGammaMass_KpKm","3#gamma mass",1600,0,4.); Pi0Gamma_kpkm= new TH1F("Pi0Gamma_kpkm","3#gamma mass",1600,0,4.); KpKmMass_vs_3GammaMass=new TH2F("KpKmMass_vs_3GammaMass","K+K- mass vs 3#gamma mass",1600,0,4.,800,0.9,2.9); PhiOmegaMass= new TH1F("PhiOmegaMass","K+K-3#gamma mass",1600,0,4.); PhiOmegaMass_with_pi0= new TH1F("PhiOmegaMass_with_pi0","K+K-3#gamma mass",1600,0,4.); } gDirectory->cd("../"); gDirectory->mkdir("KpKmP4Gamma")->cd(); { KpKm4GammaMass= new TH1F("KpKm4GammaMass","K+K-4#gamma mass",1600,0,4.); F1pi0_KpKm= new TH1F("F1pi0_KpKm","K+K-4#gamma mass",1600,0,4.); FourGammaMass_KpKm= new TH1F("FourGammaMass_KpKm","4#gamma mass",1600,0,4.); FourGamma2d_KpKm=new TH2F("FourGamma2d_KpKm","m(2#gamma,pair2) vs m(2#gamma,pair1)",100,0,1,100,0,1); TwoPi0Mass_KpKm= new TH1F("TwoPi0Mass_KpKm","4#gamma mass",1600,0,4.); TwoPi0VsKpKm= new TH2F("TwoPi0VsKpKm","2#pi0 mass vs KK mass", 400,0.9,2.9,1200,0.,3.); EtaPi0Mass_KpKm= new TH1F("EtaPi0Mass_KpKm","4#gamma mass",1600,0,4.); EtaPi0VsKpKm= new TH2F("EtaPi0VsKpKm","#eta#pi0 mass vs KK mass", 400,0.9,2.9,1200,0.5,3.5); KpKmMass_4g= new TH1F("KpKmMass_4g","K+K- mass",1600,0,4.); KpKmMass_with_2pi0= new TH1F("KpKmMass_with_2pi0","K+K- mass",1600,0,4.); KpKmMass_with_etapi0= new TH1F("KpKmMass_with_etapi0","K+K- mass",1600,0,4.); KpKmPi0_2pi0= new TH1F("KpKmPi0_2pi0","K+K-#pi0 mass",1600,0,4.); KpKmPi0_2pi0_no_phi= new TH1F("KpKmPi0_2pi0_no_phi","K+K-#pi0 mass",1600,0,4.); KpKmPi0_etapi0= new TH1F("KpKmPi0_etapi0","K+K-#pi0 mass",1600,0,4.); KstarKstar= new TH1F("KstarKstar","K+K-4#gamma mass",1600,0,4.); KpPi0_vs_KmPi0=new TH2F("KpPi0_vs_KmPi0","m(K+2#gamma,pair2) vs m(K-2#gamma,pair1)",100,0,2,100,0,2); PhiPi0Pi0Dalitz=new TH2F("PhiPi0Pi0Dalitz","M^{2}(#phi#pi^{0}_{1}) vs M^{2}(#phi#pi^{0}_{2})", 200,0,5,200,0,5); FourGamma2d_bg=new TH2F("FourGamma2d_KpKm_bg","m(2#gamma,pair2) vs m(2#gamma,pair1)",100,0,1,100,0,1); F1pi0_KpKm_bg= new TH1F("F1pi0_KpKm_bg","K+K-4#gamma mass",1600,0,4.); TwoPi0Mass_KpKm_bg= new TH1F("TwoPi0Mass_KpKm_bg","4#gamma mass",1600,0,4.); EtaPi0Mass_KpKm_bg= new TH1F("EtaPi0Mass_KpKm_bg","4#gamma mass",1600,0,4.); EtaPi0VsKpKm_bg= new TH2F("EtaPi0VsKpKm_bg","#eta#pi0 mass vs KK mass", 400,0.9,2.9,1600,0.5,3.5); KpKmMass_with_2pi0_bg= new TH1F("KpKmMass_with_2pi0_bg","K+K- mass", 800,0.9,2.9); KpKmMass_with_etapi0_bg= new TH1F("KpKmMass_with_etapi0_bg","K+K- mass", 800,0.9,2.9); KpKmPi0_2pi0_bg= new TH1F("KpKmPi0_2pi0_bg","K+K-#pi0 mass",1600,0,4.); KpKmPi0_2pi0_no_phi_bg= new TH1F("KpKmPi0_2pi0_no_phi_bg","K+K-#pi0 mass",1600,0,4.); KpKmPi0_etapi0_bg= new TH1F("KpKmPi0_etapi0_bg","K+K-#pi0 mass",1600,0,4.); KstarKstar_bg= new TH1F("KstarKstar_bg","K+K-4#gamma mass",1600,0,4.); KpPi0_vs_KmPi0_bg=new TH2F("KpPi0_vs_KmPi0_bg","m(K+2#gamma,pair2) vs m(K-2#gamma,pair1)",100,0,2,100,0,2); PhiPi0Pi0Dalitz_bg=new TH2F("PhiPi0Pi0Dalitz_bg","M^{2}(#phi#pi^{0}_{1}) vs M^{2}(#phi#pi^{0}_{2})", 200,0,5,200,0,5); PhiA0Mass= new TH1F("PhiA0Mass",";M(#Phia_{0}) [GeV];counts / 0.0025 GeV",1600,0,4.); PhiA0Mass_bg= new TH1F("PhiA0Mass_bg",";M(#Phia_{0}) [GeV];counts / 0.0025 GeV",1600,0,4.); PhiA2Mass= new TH1F("PhiA2Mass",";M(#Phia_{2}) [GeV];counts / 0.0025 GeV",1600,0,4.); PhiA2Mass_bg= new TH1F("PhiA2Mass_bg",";M(#Phia_{2}) [GeV];counts / 0.0025 GeV",1600,0,4.); } gDirectory->cd("../"); gDirectory->mkdir("KpKmP6Gamma")->cd(); { KpKm6GammaMass= new TH1F("KpKm6GammaMass","K+K-6#gamma mass",1600,0,4.); PhiEta_6g= new TH1F("PhiEta_6g","K+K-6#gamma mass",1600,0,4.); SixGammaMass_KpKm= new TH1F("SixGammaMass_KpKm","6#gamma mass",1600,0,4.); KpKmMass_6g= new TH1F("KpKmMass_6g","K+K- mass",800,0.9,2.9); KpKmEtaMass_6g_kf=new TH1F("KpKmEtaMass_6g_kf","K+K-#eta mass",1600,0,4.); KpKmEtaMass_6g_kf->SetYTitle("counts / 0.0025 GeV"); KpKmEtaMass_6g_kf->SetXTitle("M(K^{+}K^{-}#eta) [GeV]"); KpKm_vs_6g=new TH2F("KpKm_vs_6g","K+K- mass vs 6#gamma mass",600,0.,3.,200,0.9,1.4); TwoGammaMass_KpKm6g_kf= new TH1F("TwoGammaMass_KpKm6g_kf","M(2#gamma)",1600,0,4.); ThreePi0Mass_KpKm6g_kf= new TH1F("ThreePi0Mass_KpKm6g_kf","M(6#gamma)",1600,0,4.); TwoPi0EtaMass_KpKm6g_kf= new TH1F("TwoPi0EtaMass_KpKm6g_kf","M(6#gamma)",1600,0,4.); TwoEtaPi0Mass_KpKm6g_kf= new TH1F("TwoEtaPi0Mass_KpKm6g_kf","M(6#gamma)",1600,0,4.); ThreeEtaMass_KpKm6g_kf= new TH1F("ThreeEtaMass_KpKm6g_kf","M(6#gamma)",1600,0,4.); } gDirectory->cd("../"); gDirectory->mkdir("2Kp2Km")->cd(); { TwoKpTwoKmMass= new TH1F("TwoKpTwoKmMass","2K+2K- mass",1600,0,4.); PhiKpKmMass= new TH1F("PhiKpKmMass","#phiK+K- mass",800,2,4.); KpKm_vs_KpKm=new TH2F("KpKm_vs_KpKm","M(K+K- pair 1) vs M(K+K- pair 2)",800,0.9,2.9,800,0.9,2.9); P2KmMass_2kp= new TH1F("P2KmMass_2kp","p2K- mass",1600,0,4.); PKmMass_2kp= new TH1F("PKmMass_2kp","pK- mass",1600,0,4.); PKmMass_2kp_NoPhis= new TH1F("PKmMass_2kp_NoPhis","pK- mass",1600,0,4.); TwoPhiMass= new TH1F("TwoPhiMass","2K+2K- mass",1600,0,4.); TwoKpMissingMass =new TH1F("TwoKpMissingMass","M(p2K^-)",1600,1.9,4.9); TwoKpMissingMassCut =new TH1F("TwoKpMissingMassCut","M(p2K^-)",1600,1.9,4.9); } gDirectory->cd("../"); gDirectory->mkdir("KpKmPipPimP")->cd(); { PKmPipMass_vs_KpPimMass= new TH2F("PKmPipMass_vs_KpPimMass","pK-#pi+ mass vs K+#pi- mass",1600,0,4.,900,0.,4.5); KmPipMass_vs_KpPimMass= new TH2F("KmPipMass_vs_KpPimMass","K-#pi+ mass vs K+#pi- mass",1600,0,4.,900,0.,4.5); KmKpMass_vs_PipPimMass= new TH2F("KmKpMass_vs_PipPimMass","K-K+ mass vs #pi+#pi- mass",1600,0,4.,1600,0.9,4.9); KmKpMass_vs_PipPimMass->SetXTitle("M(#pi^{+}#pi^{-}) [GeV]"); KmKpMass_vs_PipPimMass->SetYTitle("M(K^{+}K^{-}) [GeV]"); PhiPipPimMass=new TH1F("PhiPipPimMass","#phi#pi+#pi- mass",1600,1.,5.); } gDirectory->cd("../"); gDirectory->mkdir("KpKmPipPim2gP")->cd(); { TwoGamma_PKpKmPipPim=new TH1F("TwoGamma_PKpKmPipPim","2#gamma mass",1600,0,4.); PKmPipMass_vs_KpPimPi0Mass= new TH2F("PKmPipMass_vs_KpPimPi0Mass","pK-#pi+ mass vs K+#pi-#pi0 mass",1600,0,4.,900,0.,4.5); KpKmMass_vs_PipPimPi0Mass= new TH2F("KpKmMass_vs_PipPimPi0Mass","K+K- mass vs #pi+#pi-#pi0 mass",1600,0,4.,800,0.9,2.9); KpKmMass_vs_PipPimPi0Mass->SetXTitle("M(#pi^{+}#pi^{-}#pi^{0}) [GeV]"); KpKmMass_vs_PipPimPi0Mass->SetYTitle("M(K^{+}K^{-}) [GeV]"); KpKmMass_vs_PipPimEtaMass= new TH2F("KpKmMass_vs_PipPimEtaMass","K+K- mass vs #pi+#pi-#eta mass",1600,0,4.,800,0.9,2.9); KpKmMass_vs_PipPimEtaMass->SetXTitle("M(#pi^{+}#pi^{-}#eta) [GeV]"); KpKmMass_vs_PipPimEtaMass->SetYTitle("M(K^{+}K^{-}) [GeV]"); PhiOmegaMass_pippim= new TH1F("PhiOmegaMass_pippim","K+K-#pi+#pi-#pi0 mass",1600,0,4.); PhiEtaMass_pippim= new TH1F("PhiEtaMass_pippim","K+K-#pi+#pi-#pi0 mass",1600,0,4.); PhiEtaPrimeMass_pippim= new TH1F("PhiEtaPrimeMass_pippim","K+K-#pi+#pi-#eta mass",1600,0,4.); } gDirectory->cd("../"); gDirectory->mkdir("KpKmPipPim4gP")->cd(); { FourGamma2d_KpKmPipPim=new TH2F("FourGamma2d_KpKmPipPim","m(2#gamma,pair2) vs m(2#gamma,pair1)",100,0,1,100,0,1); KpKmMass_with_2pi0pippim= new TH1F("KpKmMass_with_2pi0pippim","K+K- mass",800,0.9,2.9); KpKmMass_with_2pi0pippim->SetXTitle("M(K^{+}K^{-}) [GeV]"); KpKmMass_with_2pi0pippim->SetYTitle("counts / 0.0025 GeV"); PipPimPi0_with_PhiPi0=new TH1F("PipPimPi0_with_PhiPi0", ";M(#pi^{+}#pi^{-}#pi^{0}) [GeV]", 1600,0,4); PipPimPi0_with_PhiPi0_bg=new TH1F("PipPimPi0_with_PhiPi0_bg", ";M(#pi^{+}#pi^{-}#pi^{0}) [GeV]", 1600,0,4); PKmPip_vs_KpPim2Pi0=new TH2F("PKmPip_vs_KpPim2Pi0", ";M(K^{+}#pi^{-}2#pi^{0}) [GeV];M(pK^-#pi^+) [GEV]",1600,0,4.,900,0.,4.5); PKmPip_vs_KpPim2Pi0_Phicut=new TH2F("PKmPip_vs_KpPim2Pi0_Phicut", ";M(K^{+}#pi^{-}2#pi^{0}) [GeV];M(pK^-#pi^+) [GEV]",1600,0,4.,900,0.,4.5); PipPimPi0_with_PhiPi0_2pi0bg=new TH1F("PipPimPi0_with_PhiPi0_2pi0bg", ";M(#pi^{+}#pi^{-}#pi^{0}) [GeV]", 1600,0,4); PKmPip_vs_KpPim2Pi0_bg=new TH2F("PKmPip_vs_KpPim2Pi0_bg", ";M(K^{+}#pi^{-}2#pi^{0}) [GeV];M(pK^-#pi^+) [GEV]",1600,0,4.,900,0.,4.5); PKmPip_vs_KpPim2Pi0_Phicut_bg=new TH2F("PKmPip_vs_KpPim2Pi0_Phicut_bg", ";M(K^{+}#pi^{-}2#pi^{0}) [GeV];M(pK^-#pi^+) [GEV]",1600,0,4.,900,0.,4.5); KpKmMass_with_2pi0pippim_bg= new TH1F("KpKmMass_with_2pi0pippim_bg","K+K- mass",800,0.9,2.9); KpKmMass_with_2pi0pippim_bg->SetXTitle("M(K^{+}K^{-}) [GeV]"); KpKmMass_with_2pi0pippim_bg->SetYTitle("counts / 0.0025 GeV"); } gDirectory->cd("../"); gDirectory->mkdir("2Km2Pip")->cd(); { KmPipVsKmPip=new TH2F("KmPipVsKmPip","M(K-#pi+) vs M(K-#pi+)", 100,0.5,2.5,100,0.5,2.5); } gDirectory->cd("../"); gDirectory->mkdir("2KpKmPim")->cd(); { TwoKpKmPimMass= new TH1F("TwoKpKmPimMass","2K+K-pi- mass",1600,0,4.); PPimKmMass_2kp= new TH1F("PPimKmMass_2kp","p#pi-K- mass",1600,0,4.); PPimMass_2kpkm= new TH1F("PPimMass_2kpkm","p#pi- mass",400,1,3.); LambdaKmMass_2kpkm= new TH1F("LambdaKmMass_2kpkm","p#pi-K- mass",300,1,4.); } gDirectory->cd("../"); gDirectory->mkdir("2KpKmPim1G")->cd(); { TwoKpKmPimMass_1g= new TH1F("TwoKpKmPimMass_1g","2K+K-pi- mass",1600,0,4.); PPimKmMass_2kp_1g= new TH1F("PPimKmMass_2kp_1g","p#pi-K- mass",1600,0,4.); PPimMass_2kpkm_1g= new TH1F("PPimMass_2kpkm_1g","p#pi- mass",400,1,3.); LambdaKmMass_2kpkm_1g= new TH1F("LambdaKmMass_2kpkm_1g","p#pi-K- mass",300,1,4.); LambdaGamma_2kpkm_1g= new TH1F("LambdaGamma_2kpkm_1g","#Lambda#gamma mass",400,1,3.); Sigma0KmMass_2kpkm_1g= new TH1F("Sigma0KmMass_2kpkm_1g","#Sigma^{0}K- mass",300,1,4.); TwoKpKmPimMass_1g_NU= new TH1F("TwoKpKmPimMass_1g_NU","2K+K-pi- mass",1600,0,4.); PPimKmMass_2kp_1g_NU= new TH1F("PPimKmMass_2kp_1g_NU","p#pi-K- mass",1600,0,4.); PPimMass_2kpkm_1g_NU= new TH1F("PPimMass_2kpkm_1g_NU","p#pi- mass",400,1,3.); LambdaKmMass_2kpkm_1g_NU= new TH1F("LambdaKmMass_2kpkm_1g_NU","p#pi-K- mass",300,1,4.); LambdaGamma_2kpkm_1g_NU= new TH1F("LambdaGamma_2kpkm_1g_NU","#Lambda#gamma mass",400,1,3.); Sigma0KmMass_2kpkm_1g_NU= new TH1F("Sigma0KmMass_2kpkm_1g_NU","#Sigma^{0}K- mass",300,1,4.); } gDirectory->cd("../"); gDirectory->mkdir("2KpKm2PimPip")->cd(); { PipPimMass_vs_PPimMass=new TH2F("PipPimMass_vs_PPimMass", "M(#pi-#pi+) vs M(p#pi-)", 200,0.5,2.5,200,0.,2.); LambdaKmMass_2kpkm2pimpip= new TH1F("LambdaKmMass_2kpkm2pimpip","p#pi-K- mass",300,1,4.); LambdaKmMass_with_Ks= new TH1F("LambdaKmMass_with_Ks",";M(#LambdaK^{-}) [GeV]",300,1,4.); } gDirectory->cd("../"); gDirectory->mkdir("KpPim")->cd(); { PPimMass_with_K= new TH1F("PPimMass_with_K","p#pi- mass",1600,0,4.); KPimMass_with_p= new TH1F("KPimMass_with_p",";M(K^{+}#pi^{-})[GeV]",1600,0,4.); } gDirectory->cd("../"); gDirectory->mkdir("KpPimGamma")->cd(); { PPimGammaMass_with_K= new TH1F("PPimGammaMass_with_K","p#pi-#gamma mass",1600,0,4.); PPimMass_with_Kg= new TH1F("PPimMass_with_Kg","p#pi- mass",1600,0,4.); LambdaGammaMass= new TH1F("LambdaGammaMass","p#pi-#gamma mass",1600,0,4.); LambdaGammaMass_NU= new TH1F("LambdaGammaMass_NU",";M(#Lambda#gamma) [GeV]",1600,0,4.); LambdaGammaThetaVsMass= new TH2F("LambdaGammaThetaVsMass",";#Lambda#gamma mass [GeV];cos(#theta_{H})",800,0,4.,50,-1,1); LambdaGammaMassVsE= new TH2F("LambdaGammaMassVsE","p#pi-#gamma mass vs E", 48,2.8,12.4,1600,0,4.); LambdaGammaMassVst= new TH2F("LambdaGammaMassVst","p#pi-#gamma mass vs t", 2000,0.,20.,1600,0,4.); K_mom_vs_theta_1g=new TH2F("K_mom_vs_theta_1g","K+ p vs #theta",180,0,180,100,0,10); } gDirectory->cd("../"); gDirectory->mkdir("KpPim2Gamma")->cd(); { PPimGammaMass_with_Kg= new TH1F("PPimGammaMass_with_Kg","p#pi-#gamma mass",1600,0,4.); PPimMass_with_K2g= new TH1F("PPimMass_with_K2g","p#pi- mass",1600,0,4.); LambdaGammaMass_K2g= new TH1F("LambdaGammaMass_K2g","p#pi-#gamma mass",1600,0,4.); LambdaGammaMass_K2g_NU= new TH1F("LambdaGammaMass_K2g_NU","p#pi-#gamma mass",1600,0,4.); TwoGammaMass_Kppim= new TH1F("TwoGammaMass_Kppim","2#gamma mass",1600,0,4.); LambdaPi0Mass= new TH1F("LambdaPi0Mass","p#pi-#pi0 mass",1600,0,4.); LambdaPi0Mass_NU= new TH1F("LambdaPi0Mass_NU","p#pi-#pi0 mass",1600,0,4.); PPi0_with_K= new TH1F("PPi0_with_K","p#pi0 mass",1600,0,4.); LambdaEtaMass= new TH1F("LambdaEtaMass","p#pi-#eta mass",1600,0,4.); LambdaEtaMass_NU= new TH1F("LambdaEtaMass_NU","p#pi-#eta mass",1600,0,4.); Sigma0GammaMass= new TH1F("Sigma0GammaMass","p#pi-2#gamma mass",1600,0,4.); Sigma0GammaMass_NU= new TH1F("Sigma0GammaMass_NU","p#pi-2#gamma mass",1600,0,4.); SigmaPlusPiMinusMass= new TH1F("SigmaPlusPiMinusMass","p#pi-2#gamma mass",1600,0,4.); Sigma0GammaThetaVsMass= new TH2F("Sigma0GammaThetaVsMass",";#Sigma^{0}#gamma mass [GeV];cos(#theta_{H})",800,0,4.,50,-1,1); } gDirectory->cd("../"); gDirectory->mkdir("KpPim3Gamma")->cd(); { TwoGammaMass_kppim3g= new TH1F("TwoGammaMass_kppim3g","2#gamma mass",1600,0,4.); PPimMass_kppim3g= new TH1F("PPimMass_kppim3g","p#pi- mass",1600,0,4.); LambdaGammaMass_kppim3g= new TH1F("LambdaGammaMass_kppim3g","p#pi-#gamma mass",1600,0,4.); LambdaPi0Mass_kppim3g= new TH1F("LambdaPi0Mass_kppim3g","p#pi-#pi0 mass",1600,0,4.); Sigma0Pi0Mass_kppim3g= new TH1F("Sigma0Pi0Mass_kppim3g","p#pi-#pi0 mass",1600,0,4.); TwoGammaMass_kppim3g_NU= new TH1F("TwoGammaMass_kppim3g_NU","2#gamma mass",1600,0,4.); PPimMass_kppim3g_NU= new TH1F("PPimMass_kppim3g_NU","p#pi- mass",1600,0,4.); LambdaGammaMass_kppim3g_NU= new TH1F("LambdaGammaMass_kppim3g_NU","p#pi-#gamma mass",1600,0,4.); LambdaPi0Mass_kppim3g_NU= new TH1F("LambdaPi0Mass_kppim3g_NU","p#pi-#pi0 mass",1600,0,4.); Sigma0Pi0Mass_kppim3g_NU= new TH1F("Sigma0Pi0Mass_kppim3g_NU","p#pi-#pi0 mass",1600,0,4.); } gDirectory->cd("../"); gDirectory->mkdir("KpPim4Gamma")->cd(); { PPim4GammaMass_kppim4g=new TH1F("PPim4GammaMass_kppim4g","p#pi-4#gamma invariant mass",1600,1,5.); PPim4GammaMass_kppim4g->SetXTitle("M(p#pi^{-}4#gamma) [GeV]"); PPim4GammaMass_kppim4g->SetYTitle("counts / 0.0025 GeV"); PPimMass_kppim4g= new TH1F("PPimMass_kppim4g","p#pi- mass",1600,0,4.); FourGamma2d_kppim4g=new TH2F("FourGamma2d_kppim4g","m(2#gamma,pair2) vs m(2#gamma,pair1)",100,0,1,100,0,1); LambdaPi0Pi0Dalitz=new TH2F("LambdaPi0Pi0Dalitz","M^2(#Lambda#pi0_1) vs M^2(#Lambda#pi02)",300,1,16,300,1,16); LambdaEtaPi0Dalitz=new TH2F("LambdaEtaPi0Dalitz","M^2(#Lambda#eta) vs M^2(#Lambda#pi0)",300,1,16,300,1,16); FourGammaMass_kppim4g=new TH1F("FourGammaMass_kppim4g","M(4#gamma), p#pi#pi events",1600,0,4.); FourGamma2d_kppim4g->SetYTitle("M(2#gamma,pair 2) [GeV]"); FourGamma2d_kppim4g->SetXTitle("M(2#gamma,pair 1) [GeV]"); LambdaPi0Pi0Mass= new TH1F("LambdaPi0Pi0Mass","p#pi-2#pi0 mass",1600,1,5.); } gDirectory->cd("../"); gDirectory->mkdir("2Kp2Pim")->cd(); { TwoKpTwoPimMass= new TH1F("TwoKpTwoPimMass","2K+2#pi- mass",1600,0,4.); KpPim_vs_KpPim=new TH2F("KpPim_vs_KpPim","M(K+#pi- pair 1) vs M(K+#pi- pair 2)",200,0,2,200,0,2); P2PimMass_2kp= new TH1F("P2PimMass_2kp","p2#pi- mass",1600,0,4.); PPimMass_2kp= new TH1F("PPimMass_2kp","p#pi- mass",400,1,3.); LambdaPimMass= new TH1F("LambdaPimMass","p#pi-#pi- mass",400,1,3.); } gDirectory->cd("../"); gDirectory->mkdir("2Kp2Pim2g")->cd(); { TwoGammaMass_2kp2g= new TH1F("TwoGammaMass_2kp2g","2#gamma mass",1600,0,4.); PPimMass_2kp2g= new TH1F("PPimMass_2kp2g","p#pi- mass",400,1,3.); LambdaPimMass_2g= new TH1F("LambdaPimMass_2g","p#pi-#pi- mass",400,1,3.); LambdaPimPi0Mass= new TH1F("LambdaPimPi0Mass_2g","p#pi-#pi-#pi0 mass",300,1,4.); XiMinusPi0Mass= new TH1F("XiMinusPi0Mass_2g","p#pi-#pi-#pi0 mass",300,1,4.); TwoGammaMass_2kp2g_NU= new TH1F("TwoGammaMass_2kp2g_NU","2#gamma mass",1600,0,4.); PPimMass_2kp2g_NU= new TH1F("PPimMass_2kp2g_NU","p#pi- mass",400,1,3.); LambdaPimMass_2g_NU= new TH1F("LambdaPimMass_2g_NU","p#pi-#pi- mass",400,1,3.); LambdaPimPi0Mass_NU= new TH1F("LambdaPimPi0Mass_2g_NU","p#pi-#pi-#pi0 mass",300,1,4.); XiMinusPi0Mass_NU= new TH1F("XiMinusPi0Mass_2g_NU","p#pi-#pi-#pi0 mass",300,1,4.); } gDirectory->cd("../"); gDirectory->mkdir("K2Pip3Pim")->cd(); { PPim_k2pip3pim= new TH1F("PPim_k2pip3pim","; p#pi^{-} mass [GeV]",400,1,3.); PipPim2D_k2pip3pim=new TH2F("PipPim2D_k2pip3pim",";M(#pi^{+}#pi^{-}, group 1) [GEV];M(#pi^{+}#pi^{-}, group 2) [GeV]",200,0,2,200,0,2); PipPim2D_k2pip3pim_LambdaCut=new TH2F("PipPim2D_k2pip3pim_LambdaCut",";M(#pi^{+}#pi^{-}, group 1) [GEV];M(#pi^{+}#pi^{-}, group 2) [GeV]",200,0,2,200,0,2); } gDirectory->cd("../"); gDirectory->mkdir("K2PimPip")->cd(); { PPim_vs_Kpim_with_pip=new TH2F("PPim_vs_Kpim_with_pip",";M(K^{+}#pi^{-}) [GeV];M(p#pi^{-}) [GeV]",1600,0.,4.,1600,0.9,4.9); PPim_vs_pippim_with_Kp=new TH2F("PPim_vs_pippim_with_Kp",";M(#pi^{+}#pi^{-}) [GeV];M(p#pi^{-}) [GeV]",1600,0.,4.,1600,0.9,4.9); LambdaPip_vs_Kpim=new TH2F("LambdaPip_vs_Kpim",";M(K^{+}#pi^{-}) [GeV];M(#Lambda#pi^{+}) [GeV]",1600,0.,4.,1600,0.9,4.9); } gDirectory->cd("../"); gDirectory->mkdir("K2PimPip2g")->cd(); { PPim_vs_Kpim_with_pip2g=new TH2F("PPim_vs_Kpim_with_pip2g",";M(K^{+}#pi^{-}) [GeV];M(p#pi^{-}) [GeV]",1600,0.,4.,1600,0.9,4.9); PPim_vs_pippim_with_Kp2g=new TH2F("PPim_vs_pippim_with_Kp2g",";M(#pi^{+}#pi^{-}) [GeV];M(p#pi^{-}) [GeV]",1600,0.,4.,1600,0.9,4.9); LambdaPip_vs_Kpim_with_2g=new TH2F("LambdaPip_vs_Kpim_with_2g",";M(K^{+}#pi^{-}) [GeV];M(#Lambda#pi^{+}) [GeV]",1600,0.,4.,1600,0.9,4.9); TwoGamma_k2pimpip=new TH1F("TwoGamma_k2pimpip",";M(2#gamma) [GeV]",1600,0.,4.); LambdaPip_vs_KPimPi0=new TH2F("LambdaPip_vs_KpimPi0",";M(K^{+}#pi^{-}#pi^{0}) [GeV];M(#Lambda#pi^{+}) [GeV]",1600,0.,4.,1600,0.9,4.9); LambdaPipPi0_vs_KPim=new TH2F("LambdaPipPi0_vs_Kpim",";M(K^{+}#pi^{-}) [GeV];M(#Lambda#pi^{+}#pi^{0}) [GeV]",1600,0.,4.,1600,0.9,4.9); } gDirectory->cd("../"); gDirectory->mkdir("K2PimPip4g")->cd(); { PPim_vs_Kpim_with_pip4g=new TH2F("PPim_vs_Kpim_with_pip4g",";M(K^{+}#pi^{-}) [GeV];M(p#pi^{-}) [GeV]",1600,0.,4.,1600,0.9,4.9); PPim_vs_pippim_with_Kp4g=new TH2F("PPim_vs_pippim_with_Kp4g",";M(#pi^{+}#pi^{-}) [GeV];M(p#pi^{-}) [GeV]",1600,0.,4.,1600,0.9,4.9); LambdaPip_vs_Kpim_with_4g=new TH2F("LambdaPip_vs_Kpim_with_4g",";M(K^{+}#pi^{-}) [GeV];M(#Lambda#pi^{+}) [GeV]",1600,0.,4.,1600,0.9,4.9); FourGamma_k2pimpip=new TH1F("FourGamma_k2pimpip",";M(4#gamma) [GeV]",1600,0.,4.); FourGamma2d_k2pimpip=new TH2F("FourGamma2d_k2pimpip","m(2#gamma,pair2) vs m(2#gamma,pair1)",100,0,1,100,0,1); LambdaPipPi0_vs_KPimPi0=new TH2F("LambdaPipPi0_vs_KPimPi0",";M(K^{+}#pi^{-}#pi^{0}) [GeV];M(#Lambda#pi^{+}#pi^{0}) [GeV]",1600,0.,4.,1600,0.9,4.9); } gDirectory->cd("../"); gDirectory->cd("../"); MIN_BEAM_E=0.0; gPARMS->SetDefaultParameter("STRANGENESS:MIN_BEAM_E",MIN_BEAM_E); MAX_BEAM_E=12.0; gPARMS->SetDefaultParameter("STRANGENESS:MAX_BEAM_E",MAX_BEAM_E); CL_CUT=0.01; gPARMS->SetDefaultParameter("STRANGENESS:CL_CUT",CL_CUT); TRACK_CL_CUT=0.; gPARMS->SetDefaultParameter("STRANGENESS:TRACK_CL_CUT",TRACK_CL_CUT); FCAL_POS_CUT=10.0; gPARMS->SetDefaultParameter("STRANGENESS:FCAL_POS_CUT",FCAL_POS_CUT); BCAL_THRESHOLD=0.05; gPARMS->SetDefaultParameter("STRANGENESS:BCAL_THRESHOLD",BCAL_THRESHOLD); FCAL_THRESHOLD=0.1; gPARMS->SetDefaultParameter("STRANGENESS:FCAL_THRESHOLD",FCAL_THRESHOLD); UNUSED_ENERGY_CUT=0.01; gPARMS->SetDefaultParameter("STRANGENESS:UNUSED_ENERGY_CUT",UNUSED_ENERGY_CUT); // Make a container to keep track of doublets of photons MakeDoublets(); // Make a container to keep track of triplets of photons MakeTriplets(); // set up maps for keeping track of pairs of mesons vector>my_pair_vector; particle_pairs.emplace(Pi0Pi0_,my_pair_vector); particle_pairs.emplace(Pi0Eta_,my_pair_vector); particle_pairs.emplace(Pi0EtaPrime_,my_pair_vector); particle_pairs.emplace(EtaEta_,my_pair_vector); return NOERROR; } //------------------ // brun //------------------ jerror_t JEventProcessor_strangeness::brun(JEventLoop *eventLoop, int32_t runnumber) { DApplication *dapp=dynamic_cast(eventLoop->GetJApplication()); const DGeometry *geom=dapp->GetDGeometry(runnumber); double x0,y0,z0; geom->GetFCALPosition(x0,y0,z0); mFCALCenter.SetXYZ(x0,y0,z0); return NOERROR; } //------------------ // evnt //------------------ jerror_t JEventProcessor_strangeness::evnt(JEventLoop *loop, uint64_t eventnumber) { // This is called for every event. Use of common resources like writing // to a file or filling a histogram should be mutex protected. Using // loop->Get(...) to get reconstructed objects (and thereby activating the // reconstruction algorithm) should be done outside of any mutex lock // since multiple threads may call this method at the same time. // Here's an example: // // vector mydataclasses; // loop->Get(mydataclasses); // // japp->RootFillLock(this); // ... fill historgrams or trees ... // japp->RootFillUnLock(this); vectorbeamphotons; loop->Get(beamphotons); vectortagm_geom; loop->Get(tagm_geom); vectortagh_geom; loop->Get(tagh_geom); vectortracks; loop->Get(tracks); vectorneutrals; loop->Get(neutrals); vectorthrowns; loop->Get(throwns); const DMCReaction* reaction=NULL; vector triggers; loop->Get(triggers); if (throwns.size()>0){ loop->GetSingle(reaction); //loop->Get(triggers); } japp->RootWriteLock(); double mc_beam_energy=0; bool got_true_tagged=false; if (reaction!=NULL){ //printf("got here %f\n",reaction->beam.energy()); mc_beam_energy=reaction->beam.energy(); MCBeam->Fill(reaction->beam.energy()); //MCBeam->SaveAs("junk.C"); if (tagm_geom.size()&&tagh_geom.size()){ unsigned int counter=0,column=0; if ( tagm_geom[0]->E_to_column(mc_beam_energy,column)){ double Egam=0.5*(tagm_geom[0]->getEhigh(column) +tagm_geom[0]->getElow(column)); MCBeamTagged->Fill(Egam); got_true_tagged=true; //cout << Egam << " " << mc_beam_energy << endl; } else if (tagh_geom[0]->E_to_counter(mc_beam_energy,counter)){ double Egam=0.5*(tagh_geom[0]->getEhigh(counter) +tagh_geom[0]->getElow(counter)); //cout << Egam << " " << mc_beam_energy << endl; MCBeamTagged->Fill(Egam); got_true_tagged=true; } } } /* for (unsigned int i=0;ilorentzMomentum().E(); // BeamEnergy->Fill(beam_E); //cout << beamphotons[i]->lorentzMomentum().E() << endl; } */ // double t_mc=0.; vectormcprotons,mcetas,mcgammas,mcpi0s,mcKps,mcKms; if (got_true_tagged){ for (unsigned int i=0;imech==0){ switch(throwns[i]->type){ case Proton: { mcprotons.push_back(throwns[i]->lorentzMomentum()); DVector3 mom=throwns[i]->momentum(); MCProtonThetaVsP->Fill(mom.Mag(),180./M_PI*mom.Theta()); MCProtonP->Fill(mom.Mag()); double myBeam=reaction->beam.energy(); DLorentzVector beam(0,0,myBeam,myBeam); DLorentzVector target(0,0,0,ParticleMass(Proton)); DLorentzVector missing=target+beam-throwns[i]->lorentzMomentum(); MCMissingMass->Fill(missing.M()); MCMissingMassVsEbeam->Fill(myBeam,missing.M()); // printf("Thrown: p %f count %d\n",mom.Mag(),thrown_count); } break; case Eta: mcetas.push_back(throwns[i]->lorentzMomentum()); break; case Gamma: mcgammas.push_back(throwns[i]->lorentzMomentum()); break; case Pi0: mcpi0s.push_back(throwns[i]->lorentzMomentum()); break; case KPlus: mcKps.push_back(throwns[i]->lorentzMomentum()); break; case KMinus: mcKms.push_back(throwns[i]->lorentzMomentum()); break; default: break; } } } if (mcprotons.size()==1 && mcKps.size()==1 && mcKms.size()==1){ DLorentzVector KpKm=mcKps[0]+mcKms[0]; DVector3 norm=(KpKm.Vect()).Cross(reaction->beam.lorentzMomentum().Vect()); norm.SetMag(1.); DVector3 beta_lab=(1./(KpKm.E()))*KpKm.Vect(); DLorentzVector recoil=mcprotons[0]; DLorentzVector target(0,0,0,ParticleMass(Proton)); double t_mc=(recoil-target).M2(); DLorentzVector Kp=mcKps[0]; recoil.Boost(beta_lab); DVector3 zhat=recoil.Vect(); zhat.SetMag(1.); Kp.Boost(-beta_lab); DVector3 Kpdir=Kp.Vect(); Kpdir.SetMag(1.); double myBeam=reaction->beam.energy(); if (myBeam>2.8){ int index=(int)(floor(((myBeam-2.8)/0.2))); KpCosThetaH_vs_t_mc[index]->Fill(-t_mc,Kpdir.Dot(zhat)); DVector3 xhat=norm.Cross(zhat); double Kpphi=atan2(Kpdir.Dot(norm),Kpdir.Dot(xhat)); KpCosThetaVsPhi_mc[index]->Fill(Kpphi,Kpdir.Dot(zhat)); } } } if (EMULATE_TRIGGER) { loop->Get(triggers); if (triggers.size()>0){ int trig_mask=triggers[0]->Get_L1TriggerBits(); // cout << "trigger: " << trig_mask <Fill(trig_mask); if (!((trig_mask&0x1)||(trig_mask&0x4))){ japp->RootUnLock(); //cout << "Not a valid trigger: " << trig_mask <RootUnLock(); // Assign charge track by PID map>chargedParticles; FillChargedParticleVectors(tracks,chargedParticles); if (chargedParticles[Proton].size()==0) return RESOURCE_UNAVAILABLE; japp->WriteLock("mykklock"); // Find t0 (at "vertex") for event and assign list of neutral particles double t0_rf=tracks[0]->Get_BestTrackingFOM()->t0(); DVector3 vertex=tracks[0]->Get_BestTrackingFOM()->position(); // Get final state photon candidates double unused_energy=0.; vectorgammaParticles; FillGammaParticleVector(unused_energy,t0_rf,vertex,neutrals,gammaParticles); // Fill vector of beam photons to pass on the analysis along with weights // for in-time/ out-of-time vector>beamPhotons; for (unsigned int i=0;ilorentzMomentum().E()lorentzMomentum().E()>MAX_BEAM_E) continue; double dt_st=beamphotons[i]->time()-t0_rf; double weight=1.; bool got_beam_photon=false; if (fabs(dt_st)>6.012&&fabs(dt_st)<18.036){ weight=-1./6.; got_beam_photon=true; } if (fabs(dt_st)<2.004){ got_beam_photon=true; } if (got_beam_photon){ beamPhotons.push_back(make_pair(beamphotons[i],weight)); } } // Setup for performing kinematic fits //DAnalysisUtilities *dAnalysisUtilities=new DAnalysisUtilities(loop); DKinFitUtils_GlueX *dKinFitUtils = new DKinFitUtils_GlueX(loop); DKinFitter *dKinFitter = new DKinFitter(dKinFitUtils); dKinFitUtils->Reset_NewEvent(); dKinFitter->Reset_NewEvent(); for (unsigned int i=0;iGet_ConfidenceLevel(); if (CL>CL_CUT){ DLorentzVector beam_kf; map>final_kf; GetKF4vectors(dKinFitter,beam_kf,final_kf); // Particle counts unsigned int numPiPlus=final_kf[PiPlus].size(); unsigned int numPiMinus=final_kf[PiMinus].size(); unsigned int numKPlus=final_kf[KPlus].size(); unsigned int numKMinus=final_kf[KMinus].size(); unsigned int numGamma=final_kf[Gamma].size(); unsigned int numProton=final_kf[Proton].size(); unsigned int numAntiProton=final_kf[AntiProton].size(); if (numProton==1&&numAntiProton==0){ if (numKPlus==1 && numPiMinus==2 && numKMinus==0 && numPiPlus==1){ switch(numGamma){ case 0: Kp2PimPipAnalysis(unused_energy,beam_kf,final_kf,weight); break; case 2: Kp2PimPip2GammaAnalysis(unused_energy,beam_kf,final_kf,weight); break; case 4: Kp2PimPip4GammaAnalysis(unused_energy,beam_kf,final_kf,weight); break; default: break; } } if (numKPlus==1 && numPiMinus==3 && numKMinus==0 && numPiPlus==2 && numGamma==0){ Kp3Pim2PipAnalysis(unused_energy,beam_kf,final_kf,weight); } if (numKPlus==1 && numPiMinus==1 && numKMinus==0 && numPiPlus==0){ switch(numGamma){ case 0: KpPimAnalysis(unused_energy,beam_kf,final_kf,weight); break; case 1: KpPimGammaAnalysis(unused_energy,beam_kf,final_kf,weight); break; case 2: KpPim2GammaAnalysis(unused_energy,beam_kf,final_kf,weight); break; case 3: KpPim3GammaAnalysis(unused_energy,beam_kf,final_kf,weight); break; case 4: KpPim4GammaAnalysis(unused_energy,beam_kf,final_kf,weight); break; case 6: // KpPim6GammaAnalysis(unused_energy,beam_kf,final_kf,weight); break; default: break; } } if (numKPlus==2 && numPiMinus==2 && numPiPlus==0 && numKMinus==0){ switch(numGamma){ case 0: TwoKpTwoPimAnalysis(unused_energy,beam_kf,final_kf,weight); break; case 2: TwoKpTwoPim2GammaAnalysis(unused_energy,beam_kf,final_kf,weight); break; default: break; } } if (numKPlus==2&&numKMinus==2&&numPiPlus==0&&numPiMinus==0){ switch(numGamma){ case 0: TwoKpTwoKmAnalysis(unused_energy,beam_kf,final_kf,weight); break; default: break; } } if (numKPlus==1 && numKMinus==1 && numPiMinus==0 && numPiPlus==0){ switch(numGamma){ case 0: KpKmProtonAnalysis(unused_energy,beam_kf,final_kf,weight); break; case 1: KpKm1GammaProtonAnalysis(unused_energy,beam_kf,final_kf,weight); break; case 2: KpKm2GammaProtonAnalysis(unused_energy,beam_kf,final_kf,weight); break; case 3: KpKm3GammaProtonAnalysis(unused_energy,beam_kf,final_kf,weight); break; case 4: KpKm4GammaProtonAnalysis(unused_energy,beam_kf,final_kf,weight); break; case 6: KpKm6GammaProtonAnalysis(unused_energy,beam_kf,final_kf,weight); break; default: break; } } if (numKPlus==1 && numKMinus==1 && numPiPlus==1 && numPiMinus==1){ switch(numGamma){ case 0: KpKmPipPimAnalysis(unused_energy,beam_kf,final_kf,weight); break; case 2: KpKmPipPim2GammaAnalysis(unused_energy,beam_kf,final_kf,weight); break; case 4: KpKmPipPim4GammaAnalysis(unused_energy,beam_kf,final_kf,weight); break; default: break; } } if (numKPlus==2 && numKMinus==1 && numPiMinus==2 && numPiPlus==1){ switch(numGamma){ case 0: TwoKpKm2PimPipAnalysis(unused_energy,beam_kf,final_kf,weight); break; default: break; } } if (numKPlus==2 && numKMinus==1 && numPiMinus==1 && numPiPlus==0){ switch(numGamma){ case 0: TwoKpKmPimAnalysis(unused_energy,beam_kf,final_kf,weight); break; case 1: TwoKpKmPim1GammaAnalysis(unused_energy,beam_kf,final_kf,weight); break; default: break; } } if (numKPlus==0 && numKMinus==2 && numPiPlus==2 && numPiMinus==0){ if (numGamma==0){ TwoKminusTwoPiPlusAnalysis(unused_energy,beam_kf,final_kf,weight); } } } } // CL cut } } // loop over beam photons //delete dAnalysisUtilities; delete dKinFitter; delete dKinFitUtils; japp->Unlock("mykklock"); return NOERROR; } //------------------ // erun //------------------ jerror_t JEventProcessor_strangeness::erun(void) { // This is called whenever the run number changes, before it is // changed to give you a chance to clean up before processing // events from the next run number. return NOERROR; } //------------------ // fini //------------------ jerror_t JEventProcessor_strangeness::fini(void) { /* if (FILL_ROO_DATASET){ RooData->Write(); } */ // Called before program exit after event processing is finished. return NOERROR; } void JEventProcessor_strangeness::KpKmProtonAnalysis(double unused_energy,const DLorentzVector &beam_kf, map>&final_kf, double weight) const{ DLorentzVector proton_kf=final_kf[Proton][0]; DLorentzVector kp_kf=final_kf[KPlus][0]; DLorentzVector km_kf=final_kf[KMinus][0]; DLorentzVector pKm_kf=proton_kf+km_kf; DLorentzVector KpKm_kf=km_kf+kp_kf; //double t=(beam_kf-KpKm_kf).M2(); double myphi=kp_kf.Phi(); if (myphi<-M_PI) myphi+=2.*M_PI; if (myphi>M_PI) myphi-=2.*M_PI; KpThetaVsPhi_lab->Fill(180./M_PI*myphi,180./M_PI*kp_kf.Theta(),weight); KmPMass->Fill(pKm_kf.M(),weight); KpKmPDalitz->Fill(KpKm_kf.M2(),pKm_kf.M2(),weight); KpKmMass->Fill(KpKm_kf.M(),weight); KpKmMass_vs_Ebeam->Fill(beam_kf.E(),KpKm_kf.M(),weight); /* if (FILL_ROO_DATASET){ RooMass->setVal(KpKm_kf.M()); RooT->setVal(-t); RooE->setVal(beam_kf.E()); RooWeight->setVal(weight); RooData->add(RooArgSet(*RooMass,*RooT,*RooE,*RooWeight)); } */ DLorentzVector target(0,0,0,ParticleMass(Proton)); double t_mandelstam=(proton_kf-target).M2(); double t_hyperon=(pKm_kf-target).M2(); PKm_vs_t->Fill(-t_hyperon,pKm_kf.M(),weight); if (beam_kf.E()>2.8){ int index=(int)(floor(((beam_kf.E()-2.8)/0.2))); KpKm_vs_t[index]->Fill(-t_mandelstam,KpKm_kf.M(),weight); if (PHI_CUT(KpKm_kf.M())){ DVector3 norm=KpKm_kf.Vect().Cross(beam_kf.Vect()); norm.SetMag(1.); DVector3 beta_lab=(1./(KpKm_kf.E()))*KpKm_kf.Vect(); DLorentzVector recoil=proton_kf; DLorentzVector myKp=kp_kf; recoil.Boost(beta_lab); DVector3 zhat=recoil.Vect(); zhat.SetMag(1.); myKp.Boost(-beta_lab); DVector3 Kpdir=myKp.Vect(); Kpdir.SetMag(1.); KpCosThetaH_vs_t[index]->Fill(-t_mandelstam,Kpdir.Dot(zhat),weight); DVector3 xhat=norm.Cross(zhat); double Kpphi=atan2(Kpdir.Dot(norm),Kpdir.Dot(xhat)); int index2=int(-t_mandelstam/0.2); if (index2<8 && beam_kf.E()>8. && beam_kf.E()<9.){ KpCosThetaVsPhi[index2]->Fill(Kpphi,Kpdir.Dot(zhat),weight); } } } } void JEventProcessor_strangeness::KpKm1GammaProtonAnalysis(double unused_energy,const DLorentzVector &beam_kf, map>&final_kf, double weight) const{ DLorentzVector gamma_kf=final_kf[Gamma][0]; DLorentzVector kp_kf=final_kf[KPlus][0]; DLorentzVector km_kf=final_kf[KMinus][0]; DLorentzVector proton_kf=final_kf[Proton][0]; DLorentzVector pKm_kf=proton_kf+km_kf; DLorentzVector KpKm_kf=km_kf+kp_kf; DLorentzVector KpKmg_kf=KpKm_kf+gamma_kf; KpKmGammaMass->Fill(KpKmg_kf.M(),weight); KmPMass_1g->Fill(pKm_kf.M(),weight); KpKmMass_1g->Fill(KpKm_kf.M(),weight); if (PHI_CUT(KpKm_kf.M())){ PhiGammaMass->Fill(KpKmg_kf.M(),weight); } KpKmGammaMass_vs_Ebeam->Fill(beam_kf.E(),KpKmg_kf.M(),weight); } void JEventProcessor_strangeness::KpKm2GammaProtonAnalysis(double unused_energy,const DLorentzVector &beam_kf, map>&final_kf, double weight) const{ DLorentzVector gamma_kf=final_kf[Gamma][0]; DLorentzVector kp_kf=final_kf[KPlus][0]; DLorentzVector km_kf=final_kf[KMinus][0]; DLorentzVector proton_kf=final_kf[Proton][0]; DLorentzVector pKm_kf=proton_kf+km_kf; DLorentzVector KpKm_kf=km_kf+kp_kf; DLorentzVector twogam_kf=final_kf[Gamma][0]+final_kf[Gamma][1]; DLorentzVector KpKm2g_kf=KpKm_kf+twogam_kf; double twogam_mass_kf=twogam_kf.M(); DLorentzVector Kp2g=kp_kf+twogam_kf; DLorentzVector Km2g=km_kf+twogam_kf; KpKm2GammaMass->Fill(KpKm2g_kf.M(),weight); TwoGammaMass_KpKm->Fill(twogam_mass_kf,weight); KpKmMass_2g->Fill(KpKm_kf.M(),weight); KpKmMass_vs_2GammaMass->Fill(twogam_kf.M(),KpKm_kf.M(),weight); if (ETA_CUT(twogam_mass_kf,ETA_CUT_VALUE)){ KpKmEtaMass->Fill(KpKm2g_kf.M(),weight); } if (PI0_CUT(twogam_mass_kf,PI0_CUT_VALUE)){ DLorentzVector ppi0=twogam_kf+proton_kf; PPi0_vs_KpKm->Fill(KpKm_kf.M(),ppi0.M(),weight); if (ppi0.M2()>3.){ KpKmPi0Mass_baryon_cut->Fill(KpKm2g_kf.M(),weight); } KpKmPi0Mass->Fill(KpKm2g_kf.M(),weight); KpKmPi0Dalitz->Fill(Km2g.M2(),Kp2g.M2(),weight); KpKm_with_Pi0->Fill(KpKm_kf.M(),weight); if(Kp2g.M2()>0.7&&Kp2g.M2()<0.9){ KstarKm->Fill(KpKm2g_kf.M(),weight); } if(Km2g.M2()>0.7&&Km2g.M2()<0.9){ KstarKp->Fill(KpKm2g_kf.M(),weight); } } if (PHI_CUT(KpKm_kf.M())){ if (ETA_CUT(twogam_mass_kf,ETA_CUT_VALUE)){ PhiEtaMass->Fill(KpKm2g_kf.M(),weight); } if (ETAPRIME_CUT(twogam_mass_kf,ETAPRIME_CUT_VALUE)){ PhiEtaPrimeMass->Fill(KpKm2g_kf.M(),weight); } if (PI0_CUT(twogam_mass_kf,PI0_CUT_VALUE)){ PhiPi0Mass->Fill(KpKm2g_kf.M(),weight); } } } void JEventProcessor_strangeness::KpKm3GammaProtonAnalysis(double unused_energy,const DLorentzVector &beam_kf, map>&final_kf, double weight) const{ DLorentzVector kp_kf=final_kf[KPlus][0]; DLorentzVector km_kf=final_kf[KMinus][0]; DLorentzVector proton_kf=final_kf[Proton][0]; DLorentzVector pKm_kf=proton_kf+km_kf; DLorentzVector KpKm_kf=km_kf+kp_kf; DLorentzVector threegam_kf=final_kf[Gamma][0]+final_kf[Gamma][1] +final_kf[Gamma][2]; DLorentzVector KpKm3g_kf=KpKm_kf+threegam_kf; KpKm3GammaMass->Fill(KpKm3g_kf.M(),weight); ThreeGammaMass_KpKm->Fill(threegam_kf.M(),weight); DLorentzVector pair1=final_kf[Gamma][0]+final_kf[Gamma][1]; DLorentzVector pair2=final_kf[Gamma][0]+final_kf[Gamma][2]; DLorentzVector pair3=final_kf[Gamma][1]+final_kf[Gamma][2]; TwoGammaMass_kpkm3g->Fill(pair1.M(),weight); TwoGammaMass_kpkm3g->Fill(pair2.M(),weight); TwoGammaMass_kpkm3g->Fill(pair3.M(),weight); if (PI0_CUT(pair1.M(),PI0_CUT_VALUE)&&PI0_ANTI_CUT(pair2.M(),PI0_VETO_CUT) &&PI0_ANTI_CUT(pair3.M(),PI0_VETO_CUT)){ Pi0Gamma_kpkm->Fill(threegam_kf.M(),weight); if (PHI_CUT(KpKm_kf.M()) && OMEGA_CUT(threegam_kf.M())){ PhiOmegaMass_with_pi0->Fill(KpKm3g_kf.M(),weight); } } if (PI0_CUT(pair2.M(),PI0_CUT_VALUE)&&PI0_ANTI_CUT(pair1.M(),PI0_VETO_CUT) &&PI0_ANTI_CUT(pair3.M(),PI0_VETO_CUT)){ Pi0Gamma_kpkm->Fill(threegam_kf.M(),weight); if (PHI_CUT(KpKm_kf.M()) && OMEGA_CUT(threegam_kf.M())){ PhiOmegaMass_with_pi0->Fill(KpKm3g_kf.M(),weight); } } if (PI0_CUT(pair3.M(),PI0_CUT_VALUE)&&PI0_ANTI_CUT(pair2.M(),PI0_VETO_CUT) &&PI0_ANTI_CUT(pair1.M(),PI0_VETO_CUT)){ Pi0Gamma_kpkm->Fill(threegam_kf.M(),weight); if (PHI_CUT(KpKm_kf.M()) && OMEGA_CUT(threegam_kf.M())){ PhiOmegaMass_with_pi0->Fill(KpKm3g_kf.M(),weight); } } KpKmMass_vs_3GammaMass->Fill(threegam_kf.M(),KpKm_kf.M(),weight); if (PHI_CUT(KpKm_kf.M()) && threegam_kf.M()>0.73 && threegam_kf.M()<0.82){ PhiOmegaMass->Fill(KpKm3g_kf.M(),weight); } } void JEventProcessor_strangeness::KpKm4GammaProtonAnalysis(double unused_energy,const DLorentzVector &beam_kf, map>&final_kf, double weight){ DLorentzVector gamma_kf=final_kf[Gamma][0]; DLorentzVector kp_kf=final_kf[KPlus][0]; DLorentzVector km_kf=final_kf[KMinus][0]; DLorentzVector proton_kf=final_kf[Proton][0]; DLorentzVector pKm_kf=proton_kf+km_kf; DLorentzVector KpKm_kf=km_kf+kp_kf; double kpkm_mass=KpKm_kf.M(); DLorentzVector pair1=final_kf[Gamma][0]+final_kf[Gamma][1]; DLorentzVector pair2=final_kf[Gamma][2]+final_kf[Gamma][3]; DLorentzVector fourgam_kf=final_kf[Gamma][0]+final_kf[Gamma][1]+final_kf[Gamma][2]+final_kf[Gamma][3]; double fourgam_mass=fourgam_kf.M(); DLorentzVector KpKm4g_kf=KpKm_kf+fourgam_kf; double kpkm4g_mass=KpKm4g_kf.M(); KpKm4GammaMass->Fill(KpKm4g_kf.M(),weight); FourGammaMass_KpKm->Fill(fourgam_mass,weight); KpKmMass_4g->Fill(KpKm_kf.M(),weight); MakeMesonPairs(final_kf[Gamma],weight,FourGamma2d_KpKm); if (particle_pairs[Pi0Pi0_].size()>0){ TwoPi0Mass_KpKm->Fill(fourgam_mass,weight); KpKmMass_with_2pi0->Fill(kpkm_mass,weight); for (unsigned int i=0;iFill(kpkm_mass,fourgam_mass,weight); KpKmPi0_2pi0->Fill(kpkmpi0_1.M(),weight); if (kpkmpi0_1.M()>1.25 && kpkmpi0_1.M()<1.35){ F1pi0_KpKm->Fill(kpkm4g_mass,weight); } DLorentzVector pi0_2=particle_pairs[Pi0Pi0_][i].second; DLorentzVector kpkmpi0_2=KpKm_kf+pi0_2; KpKmPi0_2pi0->Fill(kpkmpi0_2.M(),weight); if (kpkmpi0_2.M()>1.25 && kpkmpi0_2.M()<1.35){ F1pi0_KpKm->Fill(kpkm4g_mass,weight); } DLorentzVector kp_pi0_1=kp_kf+pi0_1; double kp_pi0_1_mass=kp_pi0_1.M(); DLorentzVector km_pi0_2=km_kf+pi0_2; double km_pi0_2_mass=km_pi0_2.M(); KpPi0_vs_KmPi0->Fill(km_pi0_2_mass,kp_pi0_1_mass,weight); if (KSTAR_CUT(kp_pi0_1_mass)&&KSTAR_CUT(km_pi0_2_mass)){ KstarKstar->Fill(kpkm4g_mass,weight); } DLorentzVector kp_pi0_2=kp_kf+pi0_2; double kp_pi0_2_mass=kp_pi0_2.M(); DLorentzVector km_pi0_1=km_kf+pi0_1; double km_pi0_1_mass=km_pi0_1.M(); KpPi0_vs_KmPi0->Fill(km_pi0_1_mass,kp_pi0_2_mass,weight); if (KSTAR_CUT(kp_pi0_2_mass)&&KSTAR_CUT(km_pi0_1_mass)){ KstarKstar->Fill(kpkm4g_mass,weight); } if (PHI_CUT(kpkm_mass)){ PhiPi0Pi0Dalitz->Fill(kpkmpi0_1.M2(),kpkmpi0_2.M2(),weight); } else{ KpKmPi0_2pi0_no_phi->Fill(kpkmpi0_1.M(),weight); KpKmPi0_2pi0_no_phi->Fill(kpkmpi0_2.M(),weight); } } } if (particle_pairs[Pi0Eta_].size()>0){ EtaPi0Mass_KpKm->Fill(fourgam_mass,weight); KpKmMass_with_etapi0->Fill(kpkm_mass,weight); EtaPi0VsKpKm->Fill(kpkm_mass,fourgam_mass,weight); if (PHI_CUT(kpkm_mass) && fourgam_mass>0.95 && fourgam_mass<1.05){ PhiA0Mass->Fill(kpkm4g_mass,weight); } if (PHI_CUT(kpkm_mass) && fourgam_mass>1.2 && fourgam_mass<1.4){ PhiA2Mass->Fill(kpkm4g_mass,weight); } for (unsigned int i=0;iFill(kpkmpi0.M(),weight); } } /* for (unsigned int i=0;iFill(fourgam_mass,bgweight); KpKmMass_with_2pi0_bg->Fill(kpkm_mass,bgweight); DLorentzVector pi0_1=bg_particle_pairs[Pi0Pi0_][i].first; DLorentzVector kpkmpi0_1=KpKm_kf+pi0_1; KpKmPi0_2pi0->Fill(kpkmpi0_1.M(),bgweight); if (kpkmpi0_1.M()>1.25 && kpkmpi0_1.M()<1.35){ F1pi0_KpKm_bg->Fill(kpkm4g_mass,bgweight); } DLorentzVector pi0_2=bg_particle_pairs[Pi0Pi0_][i].second; DLorentzVector kpkmpi0_2=KpKm_kf+pi0_2; KpKmPi0_2pi0_bg->Fill(kpkmpi0_2.M(),bgweight); if (kpkmpi0_2.M()>1.25 && kpkmpi0_2.M()<1.35){ F1pi0_KpKm_bg->Fill(kpkm4g_mass,bgweight); } DLorentzVector kp_pi0_1=kp_kf+pi0_1; double kp_pi0_1_mass=kp_pi0_1.M(); DLorentzVector km_pi0_2=km_kf+pi0_2; double km_pi0_2_mass=km_pi0_2.M(); KpPi0_vs_KmPi0_bg->Fill(km_pi0_2_mass,kp_pi0_1_mass,bgweight); if (KSTAR_CUT(kp_pi0_1_mass)&&KSTAR_CUT(km_pi0_2_mass)){ KstarKstar_bg->Fill(kpkm4g_mass,bgweight); } DLorentzVector kp_pi0_2=kp_kf+pi0_2; double kp_pi0_2_mass=kp_pi0_2.M(); DLorentzVector km_pi0_1=km_kf+pi0_1; double km_pi0_1_mass=km_pi0_1.M(); KpPi0_vs_KmPi0_bg->Fill(km_pi0_1_mass,kp_pi0_2_mass,bgweight); if (KSTAR_CUT(kp_pi0_2_mass)&&KSTAR_CUT(km_pi0_1_mass)){ KstarKstar_bg->Fill(kpkm4g_mass,bgweight); } if (PHI_CUT(kpkm_mass)){ PhiPi0Pi0Dalitz_bg->Fill(kpkmpi0_1.M2(),kpkmpi0_2.M2(),bgweight); } else{ KpKmPi0_2pi0_no_phi_bg->Fill(kpkmpi0_1.M(),bgweight); KpKmPi0_2pi0_no_phi_bg->Fill(kpkmpi0_2.M(),bgweight); } } for (unsigned int i=0;iFill(fourgam_mass,bgweight); KpKmMass_with_etapi0_bg->Fill(kpkm_mass,bgweight); EtaPi0VsKpKm_bg->Fill(kpkm_mass,fourgam_mass,bgweight); if (PHI_CUT(kpkm_mass) && fourgam_mass>0.95 && fourgam_mass<1.05){ PhiA0Mass_bg->Fill(kpkm4g_mass,bgweight); } if (PHI_CUT(kpkm_mass) && fourgam_mass>1.2 && fourgam_mass<1.4){ PhiA2Mass_bg->Fill(kpkm4g_mass,bgweight); } DLorentzVector pi0=bg_particle_pairs[Pi0Eta_][i].first; DLorentzVector kpkmpi0=KpKm_kf+pi0; KpKmPi0_etapi0_bg->Fill(kpkmpi0.M(),bgweight); } */ } void JEventProcessor_strangeness::KpKm6GammaProtonAnalysis(double unused_energy,const DLorentzVector &beam_kf, map>&final_kf, double weight) const{ DLorentzVector kp_kf=final_kf[KPlus][0]; DLorentzVector km_kf=final_kf[KMinus][0]; DLorentzVector proton_kf=final_kf[Proton][0]; DLorentzVector pKm_kf=proton_kf+km_kf; DLorentzVector KpKm_kf=km_kf+kp_kf; DLorentzVector pair1=final_kf[Gamma][0]+final_kf[Gamma][1]; DLorentzVector pair2=final_kf[Gamma][2]+final_kf[Gamma][3]; DLorentzVector pair3=final_kf[Gamma][4]+final_kf[Gamma][5]; DLorentzVector sixgam_kf=pair1+pair2+pair3; DLorentzVector KpKm6g_kf=KpKm_kf+sixgam_kf; double sixgam_mass_kf=sixgam_kf.M(); KpKm6GammaMass->Fill(KpKm6g_kf.M(),weight); SixGammaMass_KpKm->Fill(sixgam_kf.M(),weight); KpKmMass_6g->Fill(KpKm_kf.M(),weight); KpKm_vs_6g->Fill(sixgam_kf.M(),KpKm_kf.M(),weight); if (PHI_CUT(KpKm_kf.M()) && ETA_CUT(sixgam_kf.M(),ETA_CUT_VALUE)){ PhiEta_6g->Fill(KpKm6g_kf.M(),weight); } for (unsigned int i=0;ipi0s; vectoretas; vectorunused_pair; for (unsigned int j=0;jFill(pair1.M(),weight); } if (pi0s.size()==3){ ThreePi0Mass_KpKm6g_kf->Fill(sixgam_mass_kf,weight); if (ETA_CUT(sixgam_mass_kf,ETA_CUT_VALUE)){ KpKmEtaMass_6g_kf->Fill(KpKm6g_kf.M(),weight); } } if (pi0s.size()==2 && etas.size()==1){ TwoPi0EtaMass_KpKm6g_kf->Fill(sixgam_mass_kf,weight); } if (pi0s.size()==1 && etas.size()==2){ TwoEtaPi0Mass_KpKm6g_kf->Fill(sixgam_mass_kf,weight); } if (etas.size()==3){ ThreeEtaMass_KpKm6g_kf->Fill(sixgam_mass_kf,weight); } } } void JEventProcessor_strangeness::TwoKpTwoKmAnalysis(double unused_energy,const DLorentzVector &beam_kf, map>&final_kf, double weight) const{ DLorentzVector proton_kf=final_kf[Proton][0]; DLorentzVector kp1km1=final_kf[KPlus][0]+final_kf[KMinus][0]; DLorentzVector kp2km2=final_kf[KPlus][1]+final_kf[KMinus][1]; DLorentzVector kp1km2=final_kf[KPlus][0]+final_kf[KMinus][1]; DLorentzVector kp2km1=final_kf[KPlus][1]+final_kf[KMinus][0]; DLorentzVector twokp_twokm=kp1km1+kp2km2; DLorentzVector pkm1=proton_kf+final_kf[KMinus][0]; DLorentzVector pkm2=proton_kf+final_kf[KMinus][1]; DLorentzVector pkm1km2=pkm1+final_kf[KMinus][1]; DLorentzVector twokmp=pkm1+final_kf[KMinus][1]; TwoKpMissingMass->Fill(twokmp.M(),weight); P2KmMass_2kp->Fill(pkm1km2.M(),weight); PKmMass_2kp->Fill(pkm1.M(),weight); PKmMass_2kp->Fill(pkm2.M(),weight); if ((pkm1.M()>1.5&&pkm1.M()<1.54) || (pkm2.M()>1.5&&pkm2.M()<1.54)){ TwoKpMissingMassCut->Fill(twokmp.M(),weight); } TwoKpTwoKmMass->Fill(twokp_twokm.M(),weight); KpKm_vs_KpKm->Fill(kp1km1.M(),kp2km2.M(),weight); KpKm_vs_KpKm->Fill(kp2km1.M(),kp1km2.M(),weight); if (PHI_CUT(kp1km1.M())){ PhiKpKmMass->Fill(twokp_twokm.M(),weight); } if (PHI_CUT(kp2km2.M())){ PhiKpKmMass->Fill(twokp_twokm.M(),weight); } if (PHI_CUT(kp1km2.M())){ PhiKpKmMass->Fill(twokp_twokm.M(),weight); } if (PHI_CUT(kp2km1.M())){ PhiKpKmMass->Fill(twokp_twokm.M(),weight); } if ((PHI_CUT(kp1km1.M())&&PHI_CUT(kp2km2.M())) || (PHI_CUT(kp2km1.M())&&PHI_CUT(kp1km2.M()))){ TwoPhiMass->Fill(twokp_twokm.M(),weight); } else{ PKmMass_2kp_NoPhis->Fill(pkm1.M(),weight); PKmMass_2kp_NoPhis->Fill(pkm2.M(),weight); } } void JEventProcessor_strangeness::KpKmPipPimAnalysis(double unused_energy,const DLorentzVector &beam_kf, map>&final_kf, double weight) const{ DLorentzVector kp_kf=final_kf[KPlus][0]; DLorentzVector km_kf=final_kf[KMinus][0]; DLorentzVector pip_kf=final_kf[PiPlus][0]; DLorentzVector pim_kf=final_kf[PiMinus][0]; DLorentzVector proton_kf=final_kf[Proton][0]; DLorentzVector KmPip_kf=km_kf+pip_kf; DLorentzVector pKmPip_kf=proton_kf+km_kf+pip_kf; DLorentzVector KpPim_kf=pim_kf+kp_kf; DLorentzVector KpKm_kf=kp_kf+km_kf; DLorentzVector PipPim_kf=pip_kf+pim_kf; DLorentzVector KpKmPipPim_kf=KpKm_kf+PipPim_kf; PKmPipMass_vs_KpPimMass->Fill(KpPim_kf.M(),pKmPip_kf.M(),weight); KmPipMass_vs_KpPimMass->Fill(KpPim_kf.M(),KmPip_kf.M(),weight); KmKpMass_vs_PipPimMass->Fill(PipPim_kf.M(),KpKm_kf.M(),weight); if (PHI_CUT(KpKm_kf.M())){ PhiPipPimMass->Fill(KpKmPipPim_kf.M(),weight); } } void JEventProcessor_strangeness::KpKmPipPim4GammaAnalysis(double unused_energy,const DLorentzVector &beam_kf, map>&final_kf, double weight){ DLorentzVector kp_kf=final_kf[KPlus][0]; DLorentzVector km_kf=final_kf[KMinus][0]; DLorentzVector pip_kf=final_kf[PiPlus][0]; DLorentzVector pim_kf=final_kf[PiMinus][0]; DLorentzVector proton_kf=final_kf[Proton][0]; DLorentzVector pKmPip_kf=proton_kf+km_kf+pip_kf; DLorentzVector KpPim_kf=pim_kf+kp_kf; DLorentzVector fourgam=final_kf[Gamma][0]+final_kf[Gamma][1]+final_kf[Gamma][2]+final_kf[Gamma][3]; DLorentzVector KpPim4g=KpPim_kf+fourgam; DLorentzVector KpKm_kf=km_kf+kp_kf; double kpkm_mass=KpKm_kf.M(); DLorentzVector pippim_kf=pim_kf+pip_kf; DLorentzVector pippim4g_kf=fourgam+pim_kf+pip_kf; DLorentzVector KpKmPipPim4g_kf=KpKm_kf+pippim4g_kf; MakeMesonPairs(final_kf[Gamma],weight,FourGamma2d_KpKmPipPim); for (unsigned int k=0;kFill(kpkm_mass,weight); if (PHI_CUT(kpkm_mass)){ PipPimPi0_with_PhiPi0->Fill(pippimpi0_1.M(),weight); PipPimPi0_with_PhiPi0->Fill(pippimpi0_2.M(),weight); } if (PHI_BG_CUT(kpkm_mass,NUM_SIGMA_BG)){ PipPimPi0_with_PhiPi0_bg->Fill(pippimpi0_1.M(),weight); PipPimPi0_with_PhiPi0_bg->Fill(pippimpi0_2.M(),weight); } PKmPip_vs_KpPim2Pi0->Fill(KpPim4g.M(),pKmPip_kf.M(),weight); if (kpkm_mass>1.04){ PKmPip_vs_KpPim2Pi0_Phicut->Fill(KpPim4g.M(),pKmPip_kf.M(),weight); } } /* for (unsigned int k=0;kFill(kpkm_mass,bg_weight); if (PHI_CUT(kpkm_mass)){ PipPimPi0_with_PhiPi0_2pi0bg->Fill(pippimpi0_1.M(),weight); PipPimPi0_with_PhiPi0_2pi0bg->Fill(pippimpi0_2.M(),weight); } PKmPip_vs_KpPim2Pi0_bg->Fill(KpPim4g.M(),pKmPip_kf.M(),weight); if (kpkm_mass>1.04){ PKmPip_vs_KpPim2Pi0_Phicut_bg->Fill(KpPim4g.M(),pKmPip_kf.M(),weight); } } */ } void JEventProcessor_strangeness::KpKmPipPim2GammaAnalysis(double unused_energy,const DLorentzVector &beam_kf, map>&final_kf, double weight) const{ DLorentzVector kp_kf=final_kf[KPlus][0]; DLorentzVector km_kf=final_kf[KMinus][0]; DLorentzVector pip_kf=final_kf[PiPlus][0]; DLorentzVector pim_kf=final_kf[PiMinus][0]; DLorentzVector proton_kf=final_kf[Proton][0]; DLorentzVector pKmPip_kf=proton_kf+km_kf+pip_kf; DLorentzVector KpPim_kf=pim_kf+kp_kf; DLorentzVector twogam=final_kf[Gamma][0]+final_kf[Gamma][1]; DLorentzVector KpPim2g=KpPim_kf+twogam; DLorentzVector KpKm_kf=km_kf+kp_kf; DLorentzVector pippim2g_kf=twogam+pim_kf+pip_kf; DLorentzVector KpKmPipPim2g_kf=KpKm_kf+pippim2g_kf; TwoGamma_PKpKmPipPim->Fill(twogam.M(),weight); if (ETA_CUT(twogam.M(),ETA_CUT_VALUE)){ KpKmMass_vs_PipPimEtaMass->Fill(pippim2g_kf.M(),KpKm_kf.M(),weight); if (PHI_CUT(KpKm_kf.M())&&ETAPRIME_CUT(pippim2g_kf.M(),ETAPRIME_CUT_VALUE)){ PhiEtaPrimeMass_pippim->Fill(KpKmPipPim2g_kf.M(),weight); } } if (PI0_CUT(twogam.M(),PI0_CUT_VALUE)){ PKmPipMass_vs_KpPimPi0Mass->Fill(KpPim2g.M(),pKmPip_kf.M(),weight); KpKmMass_vs_PipPimPi0Mass->Fill(pippim2g_kf.M(),KpKm_kf.M(),weight); if (PHI_CUT(KpKm_kf.M())){ if (OMEGA_CUT(pippim2g_kf.M())){ PhiOmegaMass_pippim->Fill(KpKmPipPim2g_kf.M(),weight); } if (ETA_CUT(pippim2g_kf.M(),ETA_CUT_VALUE)){ PhiEtaMass_pippim->Fill(KpKmPipPim2g_kf.M(),weight); } } } } void JEventProcessor_strangeness::TwoKpKmPimAnalysis(double unused_energy,const DLorentzVector &beam_kf, map>&final_kf, double weight) const{ DLorentzVector km_kf=final_kf[KMinus][0]; DLorentzVector pim_kf=final_kf[PiMinus][0]; DLorentzVector proton_kf=final_kf[Proton][0]; DLorentzVector ppim=proton_kf+pim_kf; DLorentzVector ppimkm=ppim+km_kf; DLorentzVector pim2kpkm=pim_kf+km_kf+final_kf[KPlus][0]+final_kf[KPlus][1]; TwoKpKmPimMass->Fill(pim2kpkm.M(),weight); PPimMass_2kpkm->Fill(ppim.M(),weight); PPimKmMass_2kp->Fill(ppimkm.M(),weight); if (LAMBDA_CUT(ppim.M())){ LambdaKmMass_2kpkm->Fill(ppimkm.M(),weight); } } void JEventProcessor_strangeness::TwoKpKmPim1GammaAnalysis(double unused_energy,const DLorentzVector &beam_kf, map>&final_kf, double weight) const{ DLorentzVector gamma_kf=final_kf[Gamma][0]; DLorentzVector pim_kf=final_kf[PiMinus][0]; DLorentzVector km_kf=final_kf[KMinus][0]; DLorentzVector proton_kf=final_kf[Proton][0]; DLorentzVector ppim=proton_kf+pim_kf; DLorentzVector ppimkm=ppim+km_kf; DLorentzVector pim2kpkm=pim_kf+km_kf+final_kf[KPlus][0]+final_kf[KPlus][1]; TwoKpKmPimMass_1g->Fill(pim2kpkm.M(),weight); PPimMass_2kpkm_1g->Fill(ppim.M(),weight); if (LAMBDA_CUT(ppim.M())){ DLorentzVector lambda_gamma=ppim+gamma_kf; LambdaGamma_2kpkm_1g->Fill(lambda_gamma.M(),weight); if (SIGMA0_CUT(lambda_gamma.M())){ DLorentzVector ppimkmg=ppimkm+gamma_kf; Sigma0KmMass_2kpkm_1g->Fill(ppimkmg.M(),weight); } LambdaKmMass_2kpkm_1g->Fill(ppimkm.M(),weight); } if (unused_energy<0.01){ TwoKpKmPimMass_1g_NU->Fill(pim2kpkm.M(),weight); PPimMass_2kpkm_1g_NU->Fill(ppim.M(),weight); if (LAMBDA_CUT(ppim.M())){ DLorentzVector lambda_gamma=ppim+gamma_kf; LambdaGamma_2kpkm_1g_NU->Fill(lambda_gamma.M(),weight); if (SIGMA0_CUT(lambda_gamma.M())){ DLorentzVector ppimkmg=ppimkm+gamma_kf; Sigma0KmMass_2kpkm_1g_NU->Fill(ppimkmg.M(),weight); } LambdaKmMass_2kpkm_1g_NU->Fill(ppimkm.M(),weight); } } } void JEventProcessor_strangeness::TwoKpKm2PimPipAnalysis(double unused_energy,const DLorentzVector &beam_kf, map>&final_kf, double weight) const{ DLorentzVector km_kf=final_kf[KMinus][0]; DLorentzVector pip_kf=final_kf[PiPlus][0]; DLorentzVector proton_kf=final_kf[Proton][0]; DLorentzVector ppim1=proton_kf+final_kf[PiMinus][0]; DLorentzVector ppim2=proton_kf+final_kf[PiMinus][1]; DLorentzVector pippim1=pip_kf+final_kf[PiMinus][0]; DLorentzVector pippim2=pip_kf+final_kf[PiMinus][1]; DLorentzVector ppim1km=ppim1+km_kf; DLorentzVector ppim2km=ppim2+km_kf; PipPimMass_vs_PPimMass->Fill(ppim1.M(),pippim2.M(),weight); PipPimMass_vs_PPimMass->Fill(ppim2.M(),pippim1.M(),weight); if (LAMBDA_CUT(ppim1.M())){ LambdaKmMass_2kpkm2pimpip->Fill(ppim1km.M(),weight); if (KSHORT_CUT(pippim2.M())){ LambdaKmMass_with_Ks->Fill(ppim1km.M(),weight); } } if (LAMBDA_CUT(ppim2.M())){ LambdaKmMass_2kpkm2pimpip->Fill(ppim2km.M(),weight); if (KSHORT_CUT(pippim1.M())){ LambdaKmMass_with_Ks->Fill(ppim2km.M(),weight); } } } void JEventProcessor_strangeness::TwoKminusTwoPiPlusAnalysis(double unused_energy,const DLorentzVector &beam_kf, map>&final_kf, double weight) const{ DLorentzVector proton_kf=final_kf[Proton][0]; DLorentzVector km1_pip1=final_kf[KMinus][0]+final_kf[PiPlus][0]; DLorentzVector km2_pip2=final_kf[KMinus][1]+final_kf[PiPlus][1]; KmPipVsKmPip->Fill(km1_pip1.M(),km2_pip2.M(),weight); DLorentzVector km1_pip2=final_kf[KMinus][0]+final_kf[PiPlus][1]; DLorentzVector km2_pip1=final_kf[KMinus][1]+final_kf[PiPlus][0]; KmPipVsKmPip->Fill(km1_pip2.M(),km2_pip1.M(),weight); } void JEventProcessor_strangeness::MakeTriplets(){ vector>pairs; pairs.push_back(make_pair(0,1)); pairs.push_back(make_pair(2,3)); pairs.push_back(make_pair(4,5)); triplets.push_back(pairs); pairs.clear(); pairs.push_back(make_pair(0,1)); pairs.push_back(make_pair(2,4)); pairs.push_back(make_pair(3,5)); triplets.push_back(pairs); pairs.clear(); pairs.push_back(make_pair(0,1)); pairs.push_back(make_pair(2,5)); pairs.push_back(make_pair(3,4)); triplets.push_back(pairs); pairs.clear(); pairs.push_back(make_pair(0,2)); pairs.push_back(make_pair(1,3)); pairs.push_back(make_pair(4,5)); triplets.push_back(pairs); pairs.clear(); pairs.push_back(make_pair(0,2)); pairs.push_back(make_pair(1,4)); pairs.push_back(make_pair(5,3)); triplets.push_back(pairs); pairs.clear(); pairs.push_back(make_pair(0,2)); pairs.push_back(make_pair(1,5)); pairs.push_back(make_pair(3,4)); triplets.push_back(pairs); pairs.clear(); pairs.push_back(make_pair(0,3)); pairs.push_back(make_pair(2,1)); pairs.push_back(make_pair(4,5)); triplets.push_back(pairs); pairs.clear(); pairs.push_back(make_pair(0,3)); pairs.push_back(make_pair(2,4)); pairs.push_back(make_pair(1,5)); triplets.push_back(pairs); pairs.clear(); pairs.push_back(make_pair(0,3)); pairs.push_back(make_pair(2,5)); pairs.push_back(make_pair(4,1)); triplets.push_back(pairs); pairs.clear(); pairs.push_back(make_pair(0,4)); pairs.push_back(make_pair(1,3)); pairs.push_back(make_pair(2,5)); triplets.push_back(pairs); pairs.clear(); pairs.push_back(make_pair(0,4)); pairs.push_back(make_pair(1,5)); pairs.push_back(make_pair(2,3)); triplets.push_back(pairs); pairs.clear(); pairs.push_back(make_pair(0,4)); pairs.push_back(make_pair(1,2)); pairs.push_back(make_pair(3,5)); triplets.push_back(pairs); pairs.clear(); pairs.push_back(make_pair(0,5)); pairs.push_back(make_pair(1,3)); pairs.push_back(make_pair(4,2)); triplets.push_back(pairs); pairs.clear(); pairs.push_back(make_pair(0,5)); pairs.push_back(make_pair(1,4)); pairs.push_back(make_pair(3,2)); triplets.push_back(pairs); pairs.clear(); pairs.push_back(make_pair(0,5)); pairs.push_back(make_pair(1,2)); pairs.push_back(make_pair(4,3)); triplets.push_back(pairs); } void JEventProcessor_strangeness::MakeDoublets(){ vector>pairs; pairs.push_back(make_pair(0,1)); pairs.push_back(make_pair(2,3)); doublets.push_back(pairs); pairs.clear(); pairs.push_back(make_pair(0,2)); pairs.push_back(make_pair(1,3)); doublets.push_back(pairs); pairs.clear(); pairs.push_back(make_pair(0,3)); pairs.push_back(make_pair(2,1)); doublets.push_back(pairs); } // Use the PIDFOM scheme to sort tracks according to particle type void JEventProcessor_strangeness::FillChargedParticleVectors(vector&tracks,map>&chargedParticles) const { Particle_t particle_list[6]={PiPlus,PiMinus,KPlus,KMinus,Proton,AntiProton}; for (unsigned int j=0;j >probabilities; for (unsigned int i=0;i<6;i++){ if ((hyp=tracks[j]->Get_Hypothesis(particle_list[i]))!=NULL){ probabilities.push_back(make_pair(hyp->Get_FOM(),particle_list[i])); } } if (probabilities.size()>0){ sort(probabilities.begin(),probabilities.end(),SortParticleProbability); Particle_t myParticle=probabilities[0].second; hyp=tracks[j]->Get_Hypothesis(myParticle); chargedParticles[myParticle].push_back(hyp->Get_TrackTimeBased()); } } } // Get list of final state photons // Create lists of photon candidates after applying some cuts on the showers. void JEventProcessor_strangeness::FillGammaParticleVector(double &unused_energy, double t0_rf, const DVector3 &vertex, vector&neutrals, vector&gammaParticles ) const{ for (unsigned int i=0;iGet_Hypothesis(Gamma); double tdiff=gamHyp->time()-t0_rf; DLorentzVector gamma_v4=gamHyp->lorentzMomentum(); ProtonGammaTimeDiff->Fill(tdiff); // Get shower info correspnding to this hypothesis const DNeutralShower *shower=gamHyp->Get_NeutralShower(); // Look for good photons bool got_photon=(fabs(tdiff)dDetectorSystem==SYS_FCAL){ if (gamma_v4.E()dQuality(shower->dBCALFCALShower))->getPosition(); pos-=mFCALCenter; if (pos.Perp()>FCAL_RADIAL_CUT) got_photon=false; if (fabs(pos.X())dDetectorSystem==SYS_BCAL){ if (gamma_v4.E()(shower->dBCALFCALShower); if (bcalshower->z>BCAL_Z_CUT) got_photon=false; } } if (got_photon){ gammaParticles.push_back(gamHyp); } else if (fabs(tdiff)>&final_kf) const{ set>myParticles=dKinFitter->Get_KinFitParticles(); set>::iterator locParticleIterator=myParticles.begin(); for(; locParticleIterator != myParticles.end(); ++locParticleIterator){ if ((*locParticleIterator)->Get_KinFitParticleType()==d_BeamParticle){ beam_kf=(*locParticleIterator)->Get_P4(); } else if ((*locParticleIterator)->Get_KinFitParticleType()==d_DetectedParticle){ Particle_t myPID=PDGtoPType((*locParticleIterator)->Get_PID()); final_kf[myPID].push_back((*locParticleIterator)->Get_P4()); } } } // Run the kinematic fitter requiring energy and momentum conservation bool JEventProcessor_strangeness::DoKinematicFit(double t0_rf,const DVector3 &vertex, const DBeamPhoton *beamphoton, map >&chargedParticles, vector&gammaParticles, DKinFitUtils_GlueX *dKinFitUtils, DKinFitter *dKinFitter) const { set> InitialParticles, FinalParticles; dKinFitter->Reset_NewFit(); shared_ptrmyBeam=dKinFitUtils->Make_BeamParticle(beamphoton); shared_ptrmyTarget=dKinFitUtils->Make_TargetParticle(Proton); InitialParticles.insert(myBeam); InitialParticles.insert(myTarget); // Vertex info DLorentzVector vert4; vert4.SetVect(vertex); vert4.SetT(t0_rf); // Charged particles Particle_t particle_list[6]={PiPlus,PiMinus,KPlus,KMinus,Proton,AntiProton}; for (unsigned int k=0;k<6;k++){ vectormyParticles=chargedParticles[particle_list[k]]; for (unsigned int m=0;mmyParticle=dKinFitUtils->Make_DetectedParticle(myParticles[m]); FinalParticles.insert(myParticle); } } // Final state photons for (unsigned int k=0;kmyGamma=dKinFitUtils->Make_DetectedParticle(22,0,0.,vert4,gammaHyp->momentum(),0.,gammaHyp->errorMatrix()); FinalParticles.insert(myGamma); } // make constraint shared_ptrlocP4Constraint = dKinFitUtils->Make_P4Constraint(InitialParticles, FinalParticles); // set constraint dKinFitter->Add_Constraint(locP4Constraint); // PERFORM THE KINEMATIC FIT return dKinFitter->Fit_Reaction(); } bool JEventProcessor_strangeness::MakeMesonPairs(vector&gammas_kf, double weight,TH2F *histo){ // Clear the current set of vectors particle_pairs[Pi0Pi0_].clear(); particle_pairs[Pi0Eta_].clear(); particle_pairs[Pi0EtaPrime_].clear(); particle_pairs[EtaEta_].clear(); // Find pairs of mesons for (unsigned int i=0;iFill(pair1_mass,pair2_mass,weight); bool eta_1_cut=ETA_CUT(pair1_mass,ETA_CUT_VALUE); bool eta_2_cut=ETA_CUT(pair2_mass,ETA_CUT_VALUE); bool pi0_1_cut=PI0_CUT(pair1_mass,PI0_CUT_VALUE); bool pi0_2_cut=PI0_CUT(pair2_mass,PI0_CUT_VALUE); if (pi0_1_cut && pi0_2_cut){ particle_pairs[Pi0Pi0_].push_back(make_pair(pair1,pair2)); } else if (eta_1_cut && eta_2_cut){ particle_pairs[EtaEta_].push_back(make_pair(pair1,pair2)); } else if (eta_1_cut && pi0_2_cut){ particle_pairs[Pi0Eta_].push_back(make_pair(pair2,pair1)); } else if (eta_2_cut && pi0_1_cut){ particle_pairs[Pi0Eta_].push_back(make_pair(pair1,pair2)); } else if (pi0_1_cut && ETAPRIME_CUT(pair2_mass,ETAPRIME_CUT_VALUE)){ particle_pairs[Pi0EtaPrime_].push_back(make_pair(pair1,pair2)); } else if (pi0_2_cut && ETAPRIME_CUT(pair1_mass,ETAPRIME_CUT_VALUE)){ particle_pairs[Pi0EtaPrime_].push_back(make_pair(pair2,pair1)); } } return (particle_pairs[Pi0Pi0_].size()>0 || particle_pairs[Pi0Eta_].size()>0 || particle_pairs[EtaEta_].size()>0 || particle_pairs[Pi0EtaPrime_].size()>0 ); } void JEventProcessor_strangeness::KpPimAnalysis(double unused_energy,const DLorentzVector &beam_kf, map>&final_kf,double weight) const { DLorentzVector ppim_kf=final_kf[Proton][0]+final_kf[PiMinus][0]; DLorentzVector kpim=final_kf[PiMinus][0]+final_kf[KPlus][0]; PPimMass_with_K->Fill(ppim_kf.M(),weight); KPimMass_with_p->Fill(kpim.M(),weight); } void JEventProcessor_strangeness::KpPimGammaAnalysis(double unused_energy,const DLorentzVector &beam_kf, map>&final_kf,double weight) const { DLorentzVector ppim_kf=final_kf[Proton][0]+final_kf[PiMinus][0]; DLorentzVector ppimgam=ppim_kf+final_kf[Gamma][0]; DLorentzVector kp_kf=final_kf[KPlus][0]; double t_mandelstam=(beam_kf-kp_kf).M2(); PPimMass_with_Kg->Fill(ppim_kf.M(),weight); PPimGammaMass_with_K->Fill(ppimgam.M(),weight); if (LAMBDA_CUT(ppim_kf.M())){ // Angular distribution in helicity frame DVector3 beta_lab=(1./(ppimgam.E()))*ppimgam.Vect(); DLorentzVector gamma_kf=final_kf[Gamma][0]; gamma_kf.Boost(-beta_lab); double cosTheta=(1./gamma_kf.E())*gamma_kf.Vect().Dot(beta_lab); LambdaGammaThetaVsMass->Fill(ppimgam.M(),cosTheta,weight); LambdaGammaMass->Fill(ppimgam.M(),weight); LambdaGammaMassVsE->Fill(beam_kf.E(),ppimgam.M(),weight); LambdaGammaMassVst->Fill(-t_mandelstam,ppimgam.M(),weight); if (SIGMA0_CUT(ppimgam.M())){ double theta=180./M_PI*kp_kf.Theta(); double p=kp_kf.P(); K_mom_vs_theta_1g->Fill(theta,p,weight); } if (unused_energy<0.01){ LambdaGammaMass_NU->Fill(ppimgam.M(),weight); } } } void JEventProcessor_strangeness::KpPim2GammaAnalysis(double unused_energy,const DLorentzVector &beam_kf, map>&final_kf,double weight) const { DLorentzVector ppim_kf=final_kf[Proton][0]+final_kf[PiMinus][0]; DLorentzVector ppimgam1=ppim_kf+final_kf[Gamma][0]; DLorentzVector ppimgam2=ppim_kf+final_kf[Gamma][1]; DLorentzVector twogam=final_kf[Gamma][0]+final_kf[Gamma][1]; DLorentzVector ppim2g=ppim_kf+twogam; DLorentzVector ppi0=final_kf[Proton][0]+twogam; PPimMass_with_K2g->Fill(ppim_kf.M(),weight); PPimGammaMass_with_Kg->Fill(ppimgam1.M(),weight); PPimGammaMass_with_Kg->Fill(ppimgam2.M(),weight); TwoGammaMass_Kppim->Fill(twogam.M(),weight); if (LAMBDA_CUT(ppim_kf.M())){ if (PI0_ANTI_CUT(twogam.M(),PI0_VETO_CUT) && ETA_ANTI_CUT(twogam.M())){ LambdaGammaMass_K2g->Fill(ppimgam1.M(),weight); if (unused_energy<0.01){ LambdaGammaMass_K2g_NU->Fill(ppimgam1.M(),weight); } DVector3 beta_lab=(1./(ppim2g.E()))*ppim2g.Vect(); if (SIGMA0_CUT(ppimgam1.M())){ Sigma0GammaMass->Fill(ppim2g.M(),weight); // Angular distribution in helicity frame final_kf[Gamma][1].Boost(-beta_lab); double cosTheta=(1./final_kf[Gamma][1].E())*final_kf[Gamma][1].Vect().Dot(beta_lab); Sigma0GammaThetaVsMass->Fill(ppim2g.M(),cosTheta,weight); if (unused_energy<0.01){ Sigma0GammaMass_NU->Fill(ppim2g.M(),weight); } } LambdaGammaMass_K2g->Fill(ppimgam2.M(),weight); if (unused_energy<0.01){ LambdaGammaMass_K2g_NU->Fill(ppimgam2.M(),weight); } if (SIGMA0_CUT(ppimgam2.M())){ Sigma0GammaMass->Fill(ppim2g.M(),weight); // Angular distribution in helicity frame final_kf[Gamma][0].Boost(-beta_lab); double cosTheta=(1./final_kf[Gamma][0].E())*final_kf[Gamma][0].Vect().Dot(beta_lab); Sigma0GammaThetaVsMass->Fill(ppim2g.M(),cosTheta,weight); if (unused_energy<0.01){ Sigma0GammaMass_NU->Fill(ppim2g.M(),weight); } } } if (PI0_CUT(twogam.M(),PI0_CUT_VALUE)){ LambdaPi0Mass->Fill(ppim2g.M(),weight); if (unused_energy<0.01){ LambdaPi0Mass_NU->Fill(ppim2g.M(),weight); } } if (ETA_CUT(twogam.M(),ETA_CUT_VALUE)){ LambdaEtaMass->Fill(ppim2g.M(),weight); if (unused_energy<0.01){ LambdaEtaMass_NU->Fill(ppim2g.M(),weight); } } } if (PI0_CUT(twogam.M(),PI0_CUT_VALUE)){ PPi0_with_K->Fill(ppi0.M(),weight); if (SIGMAP_CUT(ppi0.M())){ SigmaPlusPiMinusMass->Fill(ppim2g.M(),weight); } } } void JEventProcessor_strangeness::KpPim3GammaAnalysis(double unused_energy,const DLorentzVector &beam_kf, map>&final_kf,double weight) const { DLorentzVector ppim_kf=final_kf[Proton][0]+final_kf[PiMinus][0]; DLorentzVector ppimgam1=ppim_kf+final_kf[Gamma][0]; DLorentzVector ppimgam2=ppim_kf+final_kf[Gamma][1]; DLorentzVector ppimgam3=ppim_kf+final_kf[Gamma][2]; DLorentzVector twogam1=final_kf[Gamma][2]+final_kf[Gamma][1]; DLorentzVector twogam2=final_kf[Gamma][2]+final_kf[Gamma][0]; DLorentzVector twogam3=final_kf[Gamma][0]+final_kf[Gamma][1]; TwoGammaMass_kppim3g->Fill(twogam1.M(),weight); TwoGammaMass_kppim3g->Fill(twogam2.M(),weight); TwoGammaMass_kppim3g->Fill(twogam3.M(),weight); PPimMass_kppim3g->Fill(ppim_kf.M(),weight); if (LAMBDA_CUT(ppim_kf.M())){ LambdaGammaMass_kppim3g->Fill(ppimgam1.M(),weight); if (PI0_CUT(twogam1.M(),PI0_CUT_VALUE)){ DLorentzVector ppimpi0=twogam1+ppim_kf; LambdaPi0Mass_kppim3g->Fill(ppimpi0.M(),weight); if (SIGMA0_CUT(ppimgam1.M())){ DLorentzVector ppim3g=ppimgam1+twogam1; Sigma0Pi0Mass_kppim3g->Fill(ppim3g.M(),weight); } } LambdaGammaMass_kppim3g->Fill(ppimgam2.M(),weight); if (PI0_CUT(twogam2.M(),PI0_CUT_VALUE)){ DLorentzVector ppimpi0=twogam2+ppim_kf; LambdaPi0Mass_kppim3g->Fill(ppimpi0.M(),weight); if (SIGMA0_CUT(ppimgam2.M())){ DLorentzVector ppim3g=ppimgam2+twogam2; Sigma0Pi0Mass_kppim3g->Fill(ppim3g.M(),weight); } } LambdaGammaMass_kppim3g->Fill(ppimgam3.M(),weight); if (PI0_CUT(twogam3.M(),PI0_CUT_VALUE)){ DLorentzVector ppimpi0=twogam3+ppim_kf; LambdaPi0Mass_kppim3g->Fill(ppimpi0.M(),weight); if (SIGMA0_CUT(ppimgam3.M())){ DLorentzVector ppim3g=ppimgam3+twogam3; Sigma0Pi0Mass_kppim3g->Fill(ppim3g.M(),weight); } } } if (unused_energy<0.01){ TwoGammaMass_kppim3g_NU->Fill(twogam1.M(),weight); TwoGammaMass_kppim3g_NU->Fill(twogam2.M(),weight); TwoGammaMass_kppim3g_NU->Fill(twogam3.M(),weight); PPimMass_kppim3g_NU->Fill(ppim_kf.M(),weight); if (LAMBDA_CUT(ppim_kf.M())){ LambdaGammaMass_kppim3g_NU->Fill(ppimgam1.M(),weight); if (PI0_CUT(twogam1.M(),PI0_CUT_VALUE)){ DLorentzVector ppimpi0=twogam1+ppim_kf; LambdaPi0Mass_kppim3g_NU->Fill(ppimpi0.M(),weight); if (SIGMA0_CUT(ppimgam1.M())){ DLorentzVector ppim3g=ppimgam1+twogam1; Sigma0Pi0Mass_kppim3g_NU->Fill(ppim3g.M(),weight); } } LambdaGammaMass_kppim3g_NU->Fill(ppimgam2.M(),weight); if (PI0_CUT(twogam2.M(),PI0_CUT_VALUE)){ DLorentzVector ppimpi0=twogam2+ppim_kf; LambdaPi0Mass_kppim3g_NU->Fill(ppimpi0.M(),weight); if (SIGMA0_CUT(ppimgam2.M())){ DLorentzVector ppim3g=ppimgam2+twogam2; Sigma0Pi0Mass_kppim3g_NU->Fill(ppim3g.M(),weight); } } LambdaGammaMass_kppim3g_NU->Fill(ppimgam3.M(),weight); if (PI0_CUT(twogam3.M(),PI0_CUT_VALUE)){ DLorentzVector ppimpi0=twogam3+ppim_kf; LambdaPi0Mass_kppim3g_NU->Fill(ppimpi0.M(),weight); if (SIGMA0_CUT(ppimgam3.M())){ DLorentzVector ppim3g=ppimgam3+twogam3; Sigma0Pi0Mass_kppim3g_NU->Fill(ppim3g.M(),weight); } } } } } void JEventProcessor_strangeness::KpPim4GammaAnalysis(double unused_energy,const DLorentzVector &beam_kf, map>&final_kf,double weight){ DLorentzVector ppim_kf=final_kf[Proton][0]+final_kf[PiMinus][0]; double ppim_mass=ppim_kf.M(); DLorentzVector pair1=final_kf[Gamma][0]+final_kf[Gamma][1]; DLorentzVector pair2=final_kf[Gamma][3]+final_kf[Gamma][2]; double fourgam_mass=(pair1+pair2).M(); DLorentzVector ppim4g_kf=ppim_kf+pair1+pair2; double ppim4g_mass=ppim4g_kf.M(); PPim4GammaMass_kppim4g->Fill(ppim4g_mass,weight); FourGammaMass_kppim4g->Fill(fourgam_mass,weight); PPimMass_kppim4g->Fill(ppim_mass,weight); MakeMesonPairs(final_kf[Gamma],weight,FourGamma2d_kppim4g); if (particle_pairs[Pi0Pi0_].size()>0){ if (LAMBDA_CUT(ppim_mass)){ LambdaPi0Pi0Mass->Fill(ppim4g_mass,weight); } for (unsigned int i=0;iFill(lambdapi01.M2(),lambdapi02.M2(),weight); } } } if (particle_pairs[Pi0Eta_].size()>0){ for (unsigned int i=0;iFill(lambdapi0.M2(),lambdaeta.M2(),weight); } } } } void JEventProcessor_strangeness::TwoKpTwoPimAnalysis(double unused_energy,const DLorentzVector &beam_kf, map>&final_kf,double weight) const { DLorentzVector kp1pim1=final_kf[KPlus][0]+final_kf[PiMinus][0]; DLorentzVector kp2pim2=final_kf[KPlus][1]+final_kf[PiMinus][1]; DLorentzVector kp1pim2=final_kf[KPlus][0]+final_kf[PiMinus][1]; DLorentzVector kp2pim1=final_kf[KPlus][1]+final_kf[PiMinus][0]; DLorentzVector twokp_twopim=kp1pim1+kp2pim2; DLorentzVector ppim1=final_kf[Proton][0]+final_kf[PiMinus][0]; DLorentzVector ppim2=final_kf[Proton][0]+final_kf[PiMinus][1]; DLorentzVector ppim1pim2=ppim1+final_kf[PiMinus][1]; P2PimMass_2kp->Fill(ppim1pim2.M(),weight); PPimMass_2kp->Fill(ppim1.M(),weight); PPimMass_2kp->Fill(ppim2.M(),weight); if (LAMBDA_CUT(ppim1.M()) || LAMBDA_CUT(ppim2.M())){ LambdaPimMass->Fill(ppim1pim2.M(),weight); } TwoKpTwoPimMass->Fill(twokp_twopim.M(),weight); KpPim_vs_KpPim->Fill(kp1pim1.M(),kp2pim2.M(),weight); KpPim_vs_KpPim->Fill(kp2pim1.M(),kp1pim2.M(),weight); } void JEventProcessor_strangeness::TwoKpTwoPim2GammaAnalysis(double unused_energy,const DLorentzVector &beam_kf, map>&final_kf,double weight) const { DLorentzVector kp1pim1=final_kf[KPlus][0]+final_kf[PiMinus][0]; DLorentzVector kp2pim2=final_kf[KPlus][1]+final_kf[PiMinus][1]; DLorentzVector kp1pim2=final_kf[KPlus][0]+final_kf[PiMinus][1]; DLorentzVector kp2pim1=final_kf[KPlus][1]+final_kf[PiMinus][0]; DLorentzVector twokp_twopim=kp1pim1+kp2pim2; DLorentzVector ppim1=final_kf[Proton][0]+final_kf[PiMinus][0]; DLorentzVector ppim2=final_kf[Proton][0]+final_kf[PiMinus][1]; DLorentzVector ppim1pim2=ppim1+final_kf[PiMinus][1]; DLorentzVector twogam_kf=final_kf[Gamma][0]+final_kf[Gamma][1]; double twogam_mass_kf=twogam_kf.M(); TwoGammaMass_2kp2g->Fill(twogam_mass_kf,weight); PPimMass_2kp2g->Fill(ppim1.M(),weight); PPimMass_2kp2g->Fill(ppim2.M(),weight); if (LAMBDA_CUT(ppim1.M()) || LAMBDA_CUT(ppim2.M())){ LambdaPimMass_2g->Fill(ppim1pim2.M(),weight); if (PI0_CUT(twogam_mass_kf,PI0_CUT_VALUE)){ DLorentzVector ppim1pim2pi0=ppim1pim2+twogam_kf; LambdaPimPi0Mass->Fill(ppim1pim2pi0.M(),weight); if (XIMINUS_CUT(ppim1pim2.M())){ XiMinusPi0Mass->Fill(ppim1pim2pi0.M(),weight); } } } if (unused_energy<0.01){ TwoGammaMass_2kp2g_NU->Fill(twogam_mass_kf,weight); PPimMass_2kp2g_NU->Fill(ppim1.M(),weight); PPimMass_2kp2g_NU->Fill(ppim2.M(),weight); if (LAMBDA_CUT(ppim1.M()) || LAMBDA_CUT(ppim2.M())){ LambdaPimMass_2g_NU->Fill(ppim1pim2.M(),weight); if (PI0_CUT(twogam_mass_kf,PI0_CUT_VALUE)){ DLorentzVector ppim1pim2pi0=ppim1pim2+twogam_kf; LambdaPimPi0Mass_NU->Fill(ppim1pim2pi0.M(),weight); if (XIMINUS_CUT(ppim1pim2.M())){ XiMinusPi0Mass_NU->Fill(ppim1pim2pi0.M(),weight); } } } } } void JEventProcessor_strangeness::Kp2PimPipAnalysis(double unused_energy,const DLorentzVector &beam_kf, map>&final_kf,double weight) const{ DLorentzVector proton_kf=final_kf[Proton][0]; DLorentzVector kp_kf=final_kf[KPlus][0]; DLorentzVector pip_kf=final_kf[PiPlus][0]; DLorentzVector ppim1=proton_kf+final_kf[PiMinus][0]; DLorentzVector ppim2=proton_kf+final_kf[PiMinus][1]; DLorentzVector kpim1=kp_kf+final_kf[PiMinus][0]; DLorentzVector kpim2=kp_kf+final_kf[PiMinus][1]; DLorentzVector pippim1=pip_kf+final_kf[PiMinus][0]; DLorentzVector pippim2=pip_kf+final_kf[PiMinus][1]; PPim_vs_Kpim_with_pip->Fill(kpim1.M(),ppim2.M(),weight); PPim_vs_Kpim_with_pip->Fill(kpim2.M(),ppim1.M(),weight); PPim_vs_pippim_with_Kp->Fill(pippim2.M(),ppim1.M(),weight); PPim_vs_pippim_with_Kp->Fill(pippim1.M(),ppim2.M(),weight); if (LAMBDA_CUT(ppim1.M())){ DLorentzVector lambdapip=ppim1+pip_kf; LambdaPip_vs_Kpim->Fill(kpim2.M(),lambdapip.M(),weight); } if (LAMBDA_CUT(ppim2.M())){ DLorentzVector lambdapip=ppim2+pip_kf; LambdaPip_vs_Kpim->Fill(kpim1.M(),lambdapip.M(),weight); } } void JEventProcessor_strangeness::Kp2PimPip2GammaAnalysis(double unused_energy,const DLorentzVector &beam_kf,map>&final_kf,double weight) const{ DLorentzVector proton_kf=final_kf[Proton][0]; DLorentzVector kp_kf=final_kf[KPlus][0]; DLorentzVector pip_kf=final_kf[PiPlus][0]; DLorentzVector ppim1=proton_kf+final_kf[PiMinus][0]; DLorentzVector ppim2=proton_kf+final_kf[PiMinus][1]; DLorentzVector kpim1=kp_kf+final_kf[PiMinus][0]; DLorentzVector kpim2=kp_kf+final_kf[PiMinus][1]; DLorentzVector pippim1=pip_kf+final_kf[PiMinus][0]; DLorentzVector pippim2=pip_kf+final_kf[PiMinus][1]; DLorentzVector twogam=final_kf[Gamma][0]+final_kf[Gamma][1]; TwoGamma_k2pimpip->Fill(twogam.M(),weight); PPim_vs_Kpim_with_pip2g->Fill(kpim1.M(),ppim2.M(),weight); PPim_vs_Kpim_with_pip2g->Fill(kpim2.M(),ppim1.M(),weight); PPim_vs_pippim_with_Kp2g->Fill(pippim2.M(),ppim1.M(),weight); PPim_vs_pippim_with_Kp2g->Fill(pippim1.M(),ppim2.M(),weight); if (LAMBDA_CUT(ppim1.M())){ DLorentzVector lambdapip=ppim1+pip_kf; LambdaPip_vs_Kpim_with_2g->Fill(kpim2.M(),lambdapip.M(),weight); DLorentzVector lambdapippi0=lambdapip+twogam; DLorentzVector kpimpi0=kpim2+twogam; if (PI0_CUT(twogam.M(),PI0_CUT_VALUE)){ LambdaPip_vs_KPimPi0->Fill(kpimpi0.M(),lambdapip.M(),weight); LambdaPipPi0_vs_KPim->Fill(kpim2.M(),lambdapippi0.M(),weight); } } if (LAMBDA_CUT(ppim2.M())){ DLorentzVector lambdapip=ppim2+pip_kf; LambdaPip_vs_Kpim_with_2g->Fill(kpim1.M(),lambdapip.M(),weight); DLorentzVector lambdapippi0=lambdapip+twogam; DLorentzVector kpimpi0=kpim1+twogam; if (PI0_CUT(twogam.M(),PI0_CUT_VALUE)){ LambdaPip_vs_KPimPi0->Fill(kpimpi0.M(),lambdapip.M(),weight); LambdaPipPi0_vs_KPim->Fill(kpim1.M(),lambdapippi0.M(),weight); } } } void JEventProcessor_strangeness::Kp2PimPip4GammaAnalysis(double unused_energy,const DLorentzVector &beam_kf,map>&final_kf,double weight){ DLorentzVector proton_kf=final_kf[Proton][0]; DLorentzVector kp_kf=final_kf[KPlus][0]; DLorentzVector pip_kf=final_kf[PiPlus][0]; DLorentzVector ppim1=proton_kf+final_kf[PiMinus][0]; DLorentzVector ppim2=proton_kf+final_kf[PiMinus][1]; DLorentzVector kpim1=kp_kf+final_kf[PiMinus][0]; DLorentzVector kpim2=kp_kf+final_kf[PiMinus][1]; DLorentzVector pippim1=pip_kf+final_kf[PiMinus][0]; DLorentzVector pippim2=pip_kf+final_kf[PiMinus][1]; DLorentzVector fourgam=final_kf[Gamma][0]+final_kf[Gamma][1]+final_kf[Gamma][2]+final_kf[Gamma][3]; FourGamma_k2pimpip->Fill(fourgam.M(),weight); PPim_vs_Kpim_with_pip4g->Fill(kpim1.M(),ppim2.M(),weight); PPim_vs_Kpim_with_pip4g->Fill(kpim2.M(),ppim1.M(),weight); PPim_vs_pippim_with_Kp4g->Fill(pippim2.M(),ppim1.M(),weight); PPim_vs_pippim_with_Kp4g->Fill(pippim1.M(),ppim2.M(),weight); MakeMesonPairs(final_kf[Gamma],weight,FourGamma2d_k2pimpip); if (LAMBDA_CUT(ppim1.M())){ DLorentzVector lambdapip=ppim1+pip_kf; LambdaPip_vs_Kpim_with_4g->Fill(kpim2.M(),lambdapip.M(),weight); for (unsigned int i=0;iFill(kpimpi0_1.M(),lambdapippi0_1.M(),weight); DLorentzVector kpimpi0_2=kpim2+pi0_2; DLorentzVector lambdapippi0_2=lambdapip+pi0_1; LambdaPipPi0_vs_KPimPi0->Fill(kpimpi0_2.M(),lambdapippi0_2.M(),weight); } } if (LAMBDA_CUT(ppim2.M())){ DLorentzVector lambdapip=ppim2+pip_kf; LambdaPip_vs_Kpim_with_4g->Fill(kpim1.M(),lambdapip.M(),weight); for (unsigned int i=0;iFill(kpimpi0_1.M(),lambdapippi0_1.M(),weight); DLorentzVector kpimpi0_2=kpim1+pi0_2; DLorentzVector lambdapippi0_2=lambdapip+pi0_1; LambdaPipPi0_vs_KPimPi0->Fill(kpimpi0_2.M(),lambdapippi0_2.M(),weight); } } } void JEventProcessor_strangeness::Kp3Pim2PipAnalysis(double unused_energy,const DLorentzVector &beam_kf, map>&final_kf,double weight) const{ const DLorentzVector proton_kf=final_kf[Proton][0]; for (unsigned int i=0;i<3;i++){ DLorentzVector ppim=proton_kf+final_kf[PiMinus][i]; unsigned int ind1=0,ind2=0; switch(i){ case 0: ind1=1; ind2=2; break; case 1: ind1=0; ind2=2; break; case 2: ind1=0; ind2=1; break; default: break; } vector>>pippim_groups; vector>pippim_group; pippim_group.push_back(make_pair(final_kf[PiPlus][0],final_kf[PiMinus][ind1])); pippim_group.push_back(make_pair(final_kf[PiPlus][1],final_kf[PiMinus][ind2])); pippim_groups.push_back(pippim_group); pippim_group.clear(); pippim_group.push_back(make_pair(final_kf[PiPlus][0],final_kf[PiMinus][ind2])); pippim_group.push_back(make_pair(final_kf[PiPlus][1],final_kf[PiMinus][ind1])); pippim_groups.push_back(pippim_group); PPim_k2pip3pim->Fill(ppim.M(),weight); for (unsigned int k=0;kFill(pair1.M(),pair2.M(),weight); if (LAMBDA_CUT(ppim.M())){ PipPim2D_k2pip3pim_LambdaCut->Fill(pair1.M(),pair2.M(),weight); } } } }