Filling 2D Energy Theta Histogram

From New IAC Wiki
Jump to navigation Jump to search
 
#include <math.h>
 #include <TRandom3.h>

void Reverse3_Moller()
{
        struct evt_t
        {
                Int_t event;
                Float_t FnlMom[3], MolMom[3];
        };

        ifstream in;
        in.open("Just_Scattered_Momentum_4e8.dat");
        evt_t evt;
        Int_t nlines=0;

        TLorentzVector Fnl_Lab_4Mom;
        TLorentzVector Mol_Lab_4Mom;
        TLorentzVector Fnl_CM_4Mom;
        TLorentzVector Mol_CM_4Mom;
        TLorentzVector CMS;


        TFile *f = new TFile("Reverse3_Moller_4e8.root","RECREATE");
        TTree *tree = new TTree("Moller","Moller data from ascii file");

        TH2F* MolEThetaCM=new TH2F("MolEThetaCM","Moller E and Theta CM",42,35,56,183,90,181);

        TH1F* FnlMomLab=new TH1F("FnlMomLab","Final Scattered Electron Momentum Lab Frame ",30000,0,15000);
        TH1F* MolMomLab=new TH1F("MolMomLab","Moller Momentum Lab Frame ", 12000,0,6000);
        TH1F* FnlMomCM=new TH1F("FnlMomCM","Final Scattered Electron Momentum CM Frame ",1000,40,60);
        TH1F* MolMomCM=new TH1F("MolMomCM","Moller Momentum CM Frame ",1000,40,60);

        TH1F* FnlThetaLab=new TH1F("FinalThetaLab","Final Scattered Electron Theta Lab Frame",200,0,1);
        TH1F* MolThetaLab=new TH1F("MolThetaLab","Moller Theta Lab Frame",2000,0,90);
        TH1F* FnlThetaCM=new TH1F("FnlThetaCM","Final Scattered Electron Theta CM Frame",360,0,90);
        TH1F* MolThetaCM=new TH1F("MolThetaCM","Moller Theta CM Frame",360,90,180);

        TH1F* FnlPhiCM=new TH1F("FnlPhiCM","Final Scattered Electron Phi in CM",720,-180,180);
        TH1F* MolPhiCM=new TH1F("MolPhiCM","Final Moller Phi in CM",720,-180,180);
        TH1F* FnlPhiLab=new TH1F("FnlPhiLab","Final Scattered Electron Phi in Lab",720,-180,180);
        TH1F* MolPhiLab=new TH1F("MolPhiLab","Final Moller Phi in Lab",720,-180,180);

        tree->Branch("evt",&evt.event,"event/I:FnlPx/F:FnlPy:FnlPz:MolPx:MolPy:MolPz");

        while(in.good())
        {
                evt.event=nlines;
                in >> evt.FnlMom[0] >> evt.FnlMom[1] >> evt.FnlMom[2] >> evt.MolMom[0] >> evt.MolMom[1] >> evt.MolMom[2];
                nlines++;
                tree->Fill();

                //Define using Final Moller momentum

                        //Define 4 vector
                                Mol_Lab_4Mom.SetPxPyPzE(evt.MolMom[0],evt.MolMom[1],evt.MolMom[2],
                                        sqrt(evt.MolMom[0]*evt.MolMom[0]+evt.MolMom[1]*evt.MolMom[1]+evt.MolMom[2]*evt.MolMom[2]+.511*.511));



                //Find 4 Momentum vector componets for Final Lab frame scattered electron
                        //Define 4vector


                                Fnl_Lab_4Mom.SetPxPyPzE(evt.FnlMom[0],evt.FnlMom[1],evt.FnlMom[2],
                                        sqrt(evt.FnlMom[0]*evt.FnlMom[0]+evt.FnlMom[1]*evt.FnlMom[1]+evt.FnlMom[2]*evt.FnlMom[2]+.511*.511));
 
        
                        //Boost vectors
                                CMS=Fnl_Lab_4Mom+Mol_Lab_4Mom;
                        //Boost to CM Frame  
                                Fnl_Lab_4Mom.Boost(-CMS.BoostVector());
                                Fnl_CM_4Mom=Fnl_Lab_4Mom;
                                Mol_Lab_4Mom.Boost(-CMS.BoostVector());
                                Mol_CM_4Mom=Mol_Lab_4Mom;
        
                        //Recover the original Final 4-Momentums
                                Fnl_Lab_4Mom.SetPxPyPzE(evt.FnlMom[0],evt.FnlMom[1],evt.FnlMom[2],
                                        sqrt(evt.FnlMom[0]*evt.FnlMom[0]+evt.FnlMom[1]*evt.FnlMom[1]+evt.FnlMom[2]*evt.FnlMom[2]+.511*.511));
                                Mol_Lab_4Mom.SetPxPyPzE(evt.MolMom[0],evt.MolMom[1],evt.MolMom[2],
                                        sqrt(evt.MolMom[0]*evt.MolMom[0]+evt.MolMom[1]*evt.MolMom[1]+evt.MolMom[2]*evt.MolMom[2]+.511*.511));
        
                                //Fill 2D Histogram with Moller Energy and Theta


                                if(Mol_Lab_4Mom.E()>500)
                                {
                                        MolEThetaCM->Fill(Mol_CM_4Mom.E(),(Mol_CM_4Mom.Theta()*180/3.14159265359));
        
                                        MolThetaCM->Fill(Mol_CM_4Mom.Theta()*180/3.14159265359);
                                        FnlThetaCM->Fill(Fnl_CM_4Mom.Theta()*180/3.14159265359);
                                        MolThetaLab->Fill(Mol_Lab_4Mom.Theta()*180/3.14159265359);
                                        FnlThetaLab->Fill(Fnl_Lab_4Mom.Theta()*180/3.13159265359);   
                                        FnlPhiLab->Fill(Fnl_Lab_4Mom.Phi()*180/3.14159265359);
                                        MolPhiLab->Fill(Mol_Lab_4Mom.Phi()*180/3.14159265359);
        
        
                                        MolMomLab->Fill(Mol_Lab_4Mom.P());
                                        FnlMomLab->Fill(Fnl_Lab_4Mom.P());
                                        MolMomCM->Fill(Mol_CM_4Mom.P());
                                        FnlMomCM->Fill(Fnl_CM_4Mom.P());
                                        FnlPhiCM->Fill(Fnl_CM_4Mom.Phi()*180/3.14159265359);
                                        MolPhiCM->Fill(Mol_CM_4Mom.Phi()*180/3.14159265359);
                                }
                }//End input loop
        
        in.close();
        tree->Print();  
        tree->Write(); 
        f->Write();
        delete tree;
        delete f;
}//End main



Links

Back