]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - EMCAL/AliEMCALJetFinderPlots.cxx
Corrected posting and cleaning of digitizers (T.Kuhr)
[u/mrichter/AliRoot.git] / EMCAL / AliEMCALJetFinderPlots.cxx
index aaf614fad92a497a03d464b42ecc23e1b21b8ee8..d57d7b4d60bfe2ba0c0742686f112a629edb34bd 100755 (executable)
@@ -77,17 +77,34 @@ fhInputOutput=0;
 //     TH2F                            *fhInputOutput;  //("hJetEtRatio2","Ratio of Second Highest to Highest",100,0,1);
      
 fhRecoBinPt=0;        // ("fhRecoBinPt","Reconstructed Pt Distribution",100,0,1);
+fhRecoBinPtNoBg=0;            // ("fhRecoBinPt","Reconstructed Pt Distribution",100,0,1);
 fhRecoBinPartonPt=0;    // ("fhRecoBinPartonPt","Input Pt Distribution",100,0,1);
 fhRecoBinJetEt=0;       // ("fhRecoJetEt","E_{T}^{reco}",250,0.,250.);
 fhRecoBinInputJetEt=0;  // ("fhRecoInputJetEt","E_{T}^{reco}",250,0.,250.);
+fhJetPT =0;// new TH1F("hJetPT","P_{T} Distribution",200,0,200);
+fhPartonPT =0;// new TH1F("hPartonPT","Parton P_{T} Distribution",200,0,1);
+fhJetPT2 =0;// new TH1F("hJetPT","P_{T} Distribution",200,0,200);
+fhPartonPT2 =0;// new TH1F("hPartonPT","Parton P_{T} Distribution",200,0,1);
+fhRecoBinFragmFcn =0;//new TH1F("fhRecoBinFragmFcn","Reconstructed Frag. Fcn",100,0,1);
+fhRecoBinFragmFcnNoBg =0;//new TH1F("fhRecoBinFragmFcn","Reconstructed Frag. Fcn",100,0,1);
+fhRecoBinPartonFragmFcn =0;// new TH1F("fhRecoBinPartonFragmFcn","Input Bin Fragm Fcn Distribution",100,0,1);
+fhJetInvE=0;// = new TH1F("fhJetInvE","#frac{1}{E_{R}}",100,0,1);
+fhJetInvE2=0;// = new TH1F("fhJetInvE","#frac{1}{E_{R}}",100,0,1);
+fScaleFactor = 1.0/0.6731;
+fhBackHisto=0;
+
 }
 
 void AliEMCALJetFinderPlots::InitPlots()
 {
 //========================= CASE 1 =======================================     
-    fhFragmFcn = new TH1F("hFragmFcn","Fragmentation Function",100,0,1);
+    fhFragmFcn = new TH1F("hFragmFcn","Fragmentation Function",200,0,2);
     fhFragmFcn->Sumw2();
-    fhPartonFragmFcn = new TH1F("hPartonFragmFcn","Parton Fragmentation Function",100,0,1);
+    fhJetPT = new TH1F("hJetPT","P_{T} Distribution",200,0,200);
+    fhJetPT->Sumw2();
+    fhPartonPT = new TH1F("hPartonPT","Parton P_{T} Distribution",200,0,200);
+    fhPartonPT->Sumw2();
+    fhPartonFragmFcn = new TH1F("hPartonFragmFcn","Parton Fragmentation Function",200,0,2);
     fhPartonFragmFcn->Sumw2();
     fhPartonJT = new TH1F("hPartonJT","Track Momentum Perpendicular to Parton Axis",100,0.,10.);
     fhPartonJT->Sumw2();
@@ -99,6 +116,8 @@ void AliEMCALJetFinderPlots::InitPlots()
     fhJetPL->Sumw2();
     fhJetEt = new TH1F("hJetEt","E_{T}^{reco}",250,0.,250.);
     fhJetEt->Sumw2();
+    fhJetEtDiff = new TH1F("hJetEtDiff","E_{T}^{reco}-E_{T}^{Parton}",250,-124.,125.);
+    fhJetEtDiff->Sumw2();
     fhJetEta = new     TH1F("hJetEta","#eta_{jet}^{reco}",180,-0.9,0.9);
     fhJetEta->Sumw2();
     fhJetPhi = new     TH1F("hJetPhi","#phi_{jet}^{reco}",62,0.,3.1);
@@ -154,9 +173,13 @@ void AliEMCALJetFinderPlots::InitPlots()
 //======================= CASE 2 ======================================
   
 
-fhFragmFcn2            = new TH1F("hFragmFcn2","Fragmentation Function",100,0,1);
+fhFragmFcn2            = new TH1F("hFragmFcn2","Fragmentation Function",200,0,2);
 fhFragmFcn2->Sumw2();
-fhPartonFragmFcn2      = new TH1F("hPartonFragmFcn2","Parton Fragmentation Function",100,0,1);
+    fhJetPT2 = new TH1F("hJetPT2","P_{T} Distribution",200,0,200);
+    fhJetPT2->Sumw2();
+    fhPartonPT2 = new TH1F("hPartonPT2","Parton P_{T} Distribution",200,0,1);
+    fhPartonPT2->Sumw2();
+fhPartonFragmFcn2      = new TH1F("hPartonFragmFcn2","Parton Fragmentation Function",200,0,2);
 fhPartonFragmFcn2->Sumw2();
 fhPartonJT2            = new TH1F("hPartonJT2","Track Momentum Perpendicular to Parton Axis",100,0.,10.);
 fhPartonJT2->Sumw2();
@@ -168,6 +191,8 @@ fhJetPL2                    = new TH1F("hJetPL2","Track Momentum Parallel to Jet Axis ",100,0.,10
 fhJetPL2->Sumw2();
 fhJetEt2                       = new TH1F("hJetEt2","E_{T}^{reco}",250,0.,250.);
 fhJetEt2->Sumw2();
+    fhJetEtDiff2 = new TH1F("hJetEtDiff2","E_{T}^{reco}-E_{T}^{Parton}",250,-124.,125.);
+    fhJetEtDiff2->Sumw2();
 fhJetEta2              = new TH1F("hJetEta2","#eta_{jet}^{reco}",180,-0.9,0.9);
 fhJetEta2->Sumw2();
 fhJetPhi2              = new TH1F("hJetPhi2","#phi_{jet}^{reco}",62,0.,3.1);
@@ -195,15 +220,31 @@ fhInputOutput=    new TH2F("hInputOutput","Input and Reconstruction Correlations;In
 
 //============================== Reconstruction Bin Comparison  ============================================
 
-fhRecoBinPt =new TH1F("fhRecoBinPt","Reconstructed Pt Distribution",100,0,1);
+fhRecoBinPt =new TH1F("fhRecoBinPt","Reconstructed Pt Distribution",200,0,200);
 fhRecoBinPt->Sumw2();
-fhRecoBinPartonPt = new TH1F("fhRecoBinPartonPt","Input Pt Distribution",100,0,1);
+fhRecoBinPtNoBg =new TH1F("fhRecoBinPtNoBg","Reconstructed Pt Distribution Background Subtracted",200,0,200);
+fhRecoBinPtNoBg->Sumw2();
+fhRecoBinPartonPt = new TH1F("fhRecoBinPartonPt","Input Pt Distribution",200,0,200);
 fhRecoBinPartonPt->Sumw2();
+fhRecoBinFragmFcn =new TH1F("fhRecoBinFragmFcn","Reconstructed Frag. Fcn",200,0,2);
+fhRecoBinFragmFcn->Sumw2();
+fhRecoBinFragmFcnNoBg =new TH1F("fhRecoBinFragmFcnNoBg","Reconstructed Frag. Fcn With Background Removed",200,0,2);
+fhRecoBinFragmFcnNoBg->Sumw2();
+fhRecoBinPartonFragmFcn = new TH1F("fhRecoBinPartonFragmFcn","Input Bin Fragm Fcn Distribution",200,0,2);
+fhRecoBinPartonFragmFcn->Sumw2();
 fhRecoBinJetEt = new TH1F("fhRecoJetEt","E_{T}^{reco}",250,0.,250.);
 fhRecoBinJetEt->Sumw2();
 fhRecoBinInputJetEt = new TH1F("fhRecoInputJetEt","E_{T}^{reco}",250,0.,250.);
 fhRecoBinInputJetEt->Sumw2();
-                                
+        
+
+fhJetInvE = new TH1F("fhJetInvE","#frac{1}{E_{R}}",100,0,1);
+fhJetInvE->Sumw2();
+fhJetInvE2 = new TH1F("fhJetInvE2","#frac{1}{E_{R}}",100,0,1);
+fhJetInvE2->Sumw2();
+               
+
+
   fInitialised = kTRUE;        
   
 }
@@ -218,6 +259,7 @@ delete    fhPartonPL;// = new TH1F("hPartonPL","Track Momentum Parallel to Parto
 delete    fhJetJT;// = new TH1F("hJetJT","Track Momentum Perpendicular to Jet Axis",100,0.,10.);
 delete    fhJetPL;// = new TH1F("hJetPL","Track Momentum Parallel to Jet Axis ",100,0.,100.);
 delete    fhJetEt;// = new TH1F("hJetEt","E_{T}^{reco}",250,0.,250.);
+delete fhJetEtDiff;    // ("hJetEt2","E_{T}^{reco}",250,0.,250.);
 delete    fhJetEta;// = new       TH1F("hJetEta","#eta_{jet}^{reco}",180,-0.9,0.9);
 delete    fhJetPhi;// = new       TH1F("hJetPhi","#phi_{jet}^{reco}",62,0.,3.1);
 delete    fhPartonEta;// = new    TH1F("hPartonEta","#eta_{Parton}",180,-0.9,0.9);
@@ -234,6 +276,7 @@ delete        fhEtaPhiSpread;
        delete                          fhJetJT2;       // ("hJetJT2","Track Momentum Perpendicular to Jet Axis",100,0.,10.);
        delete                          fhJetPL2;       // ("hJetPL2","Track Momentum Parallel to Jet Axis ",100,0.,100.);
        delete                          fhJetEt2;       // ("hJetEt2","E_{T}^{reco}",250,0.,250.);
+       delete                          fhJetEtDiff2;   // ("hJetEt2","E_{T}^{reco}",250,0.,250.);
        delete                          fhJetEta2;      // ("hJetEta2","#eta_{jet}^{reco}",180,-0.9,0.9);
        delete                          fhJetPhi2;      // ("hJetPhi2","#phi_{jet}^{reco}",62,0.,3.1);
        delete                          fhPartonEta2;   // ("hPartonEta2","#eta_{Parton}",180,-0.9,0.9);
@@ -248,19 +291,38 @@ delete      fhEtaPhiSpread;
        delete                          fhEtaPhiDist2;  //("hEtaPhiDist2","Angular Distance Between First and Second",100,0,3);
 
        delete fhRecoBinPt;          // ("fhRecoBinPt","Reconstructed Pt Distribution",100,0,1);
+       delete fhRecoBinPtNoBg;          // ("fhRecoBinPt","Reconstructed Pt Distribution",100,0,1);
        delete fhRecoBinPartonPt;    // ("fhRecoBinPartonPt","Input Pt Distribution",100,0,1);
        delete fhRecoBinJetEt;       // ("fhRecoJetEt","E_{T}^{reco}",250,0.,250.);
        delete fhRecoBinInputJetEt;  // ("fhRecoInputJetEt","E_{T}^{reco}",250,0.,250.);
                                        
+       delete fhJetPT ;// new TH1F("hJetPT","P_{T} Distribution",200,0,200);
+       delete fhPartonPT ;// new TH1F("hPartonPT","Parton P_{T} Distribution",200,0,1);
+       delete fhJetPT2 ;// new TH1F("hJetPT","P_{T} Distribution",200,0,200);
+       delete fhPartonPT2 ;// new TH1F("hPartonPT","Parton P_{T} Distribution",200,0,1);
+       delete fhRecoBinFragmFcn;//new TH1F("fhRecoBinFragmFcn","Reconstructed Frag. Fcn",100,0,1);
+       delete fhRecoBinFragmFcnNoBg;//new TH1F("fhRecoBinFragmFcn","Reconstructed Frag. Fcn",100,0,1);
+       delete fhRecoBinPartonFragmFcn;// new TH1F("fhRecoBinPartonFragmFcn","Input Bin Fragm Fcn Distribution",100,0,1);
+    delete fhJetInvE;// = new TH1F("fhJetInvE","#frac{1}{E_{R}}",100,0,1);
+    delete fhJetInvE2;// = new TH1F("fhJetInvE","#frac{1}{E_{R}}",100,0,1);
 
 }      
 
-void AliEMCALJetFinderPlots::FillFromOutput(AliEMCALJetFinderOutput* output)
+void AliEMCALJetFinderPlots::FillFromOutput(AliEMCALJetFinderOutput* output, Float_t weight)
 {
        // Fill histograms from an output object
 if (!fInitialised) InitPlots();        
   fOutput = output;
   if (!fOutput) return;
+// Make some temporary histograms to make sure we subtract 
+  // background properly
+/*
+tempFragmFcnNoBg =new TH1F("tempFragmFcnNoBg","Reconstructed Frag. Fcn With Background Removed",200,0,2);
+tempPtNoBg =new TH1F("tempPtNoBg","Reconstructed Frag. Fcn With Background Removed",200,0,200);
+tempFragmFcnNoBg->Fill(count/(jethighest->Energy()*fScaleFactor),-fhBackHisto->GetBinContent(count));
+tempPtNoBg->AddBinContent(count,-fhBackHisto->GetBinContent(count));
+*/
+  
   fhNJets->Fill(fOutput->GetNJets());
   Bool_t doesJetMeetBinCriteria = 0;
   AliEMCALJet* jethighest=0; 
@@ -290,10 +352,11 @@ if (!fInitialised) InitPlots();
   }
        
   
+Float_t et=0;
 if (numappjet >=1)
 {
        Float_t theta = 2.0*atan(exp(-fOutput->GetParton(0)->Eta()));
-       Float_t et    = fOutput->GetParton(0)->Energy() * TMath::Sin(theta);
+       et    = fOutput->GetParton(0)->Energy() * TMath::Sin(theta);
        if (fOutput->GetNJets()>1)
        {
                for (Int_t counter = 0; counter<numappjet;counter++)  
@@ -333,11 +396,11 @@ if (numappjet >=1)
                }
                                                
        }
-       if ( 95.0 < jethighest->Energy() && jethighest->Energy() < 105.0  )
+       if ( 95.0 < jethighest->Energy()*fScaleFactor && jethighest->Energy()*fScaleFactor < 105.0  )
        {
                doesJetMeetBinCriteria = 1;
-               fhRecoBinJetEt->Fill(jethighest->Energy());
-               fhRecoBinInputJetEt->Fill(et);
+               fhRecoBinJetEt->Fill(jethighest->Energy()*fScaleFactor,weight);
+               fhRecoBinInputJetEt->Fill(et,weight);
        }
        fhInputOutput->Fill(et,jethighest->Energy());
                
@@ -351,23 +414,25 @@ if (numappjet > 1)
   AliEMCALParton* parton;
 
   // End finding highest and second highest and continue
-  fhJetEt2->Fill(jethighest->Energy());
-  fhJetEta2->Fill(jethighest->Eta()  );
-  fhJetPhi2->Fill(jethighest->Phi() );
+  fhJetEt2->Fill(jethighest->Energy()*fScaleFactor,weight);
+  fhJetEtDiff2->Fill(jethighest->Energy()*fScaleFactor-et,weight);
+  fhJetInvE2->Fill(1.0/(jethighest->Energy()*fScaleFactor),weight);
+  fhJetEta2->Fill(jethighest->Eta(),weight  );
+  fhJetPhi2->Fill(jethighest->Phi(),weight );
   if (nPartons ==0) return;
   parton = fOutput->GetParton(0);
   
-  fhPartonEta2->Fill( parton->Eta() );
-  fhPartonPhi2->Fill( parton->Phi() );
+  fhPartonEta2->Fill( parton->Eta(),weight );
+  fhPartonPhi2->Fill( parton->Phi(),weight );
 
   //hJetEtDiff->Fill( jet->Energy() - parton->Energy() );
-  fhEtaDiff2->Fill( jethighest->Eta() - parton->Eta() );
-  fhPhiDiff2->Fill( jethighest->Phi() - parton->Phi() );
+  fhEtaDiff2->Fill( jethighest->Eta() - parton->Eta(),weight );
+  fhPhiDiff2->Fill( jethighest->Phi() - parton->Phi(),weight );
   fhEtaPhiSpread2->Fill(jethighest->Eta()-parton->Eta(),jethighest->Phi() - parton->Phi());
-  fhJetEtSecond2->Fill(jetsecond->Energy()); 
-  fhJetEtRatio2->Fill(jetsecond->Energy()/jethighest->Energy()); 
+  fhJetEtSecond2->Fill(jetsecond->Energy()*fScaleFactor,weight); 
+  fhJetEtRatio2->Fill(jetsecond->Energy()/jethighest->Energy(),weight); 
   fhEtaPhiDist2->Fill( TMath::Sqrt((jethighest->Eta() - jetsecond->Eta())*(jethighest->Eta() - jetsecond->Eta())
-                     + (jethighest->Phi() - jetsecond->Phi())*(jethighest->Phi() - jetsecond->Phi())     ));  
+                     + (jethighest->Phi() - jetsecond->Phi())*(jethighest->Phi() - jetsecond->Phi())     ),weight);  
   /* 
   Float_t *pt,*phi,*eta;
   Int_t *pdg;
@@ -408,19 +473,20 @@ if (numappjet > 1)
                  alpha = TMath::ACos(crp*ctp*srt*stt+srp*stp*srt*stt+crt*ctt);   
          }
       Double_t correctp = pt[iT]/stt;  
-      fhPartonPL2->Fill( correctp*cos(alpha));
+      fhPartonPL2->Fill( correctp*cos(alpha),weight);
       if ( (parton->Eta()-eta[iT])*(parton->Eta()-eta[iT]) +  
       (parton->Phi()-phi[iT])*(parton->Phi()-phi[iT]) < 0.2*0.2 ) 
-             fhPartonJT2->Fill( correctp*sin(alpha));
+             fhPartonJT2->Fill( correctp*sin(alpha),weight);
+      fhPartonPT2->Fill(correctp*sin(tt),weight);
       if (fNominalEnergy == 0.0) {
-             fhPartonFragmFcn2->Fill(  correctp*sin(tt)/parton->Energy()  );
+             fhPartonFragmFcn2->Fill(  correctp*sin(tt)/parton->Energy(),weight  );
       }else 
       {
-              fhPartonFragmFcn2->Fill(correctp*sin(tt)/fNominalEnergy);
+              fhPartonFragmFcn2->Fill(correctp*sin(tt)/fNominalEnergy,weight);
       }          
       if (doesJetMeetBinCriteria)
       {
-             fhRecoBinPartonPt->Fill(correctp*sin(tt));          // ("fhRecoBinPt","Reconstructed Pt Distribution",100,0,1);
+             fhRecoBinPartonPt->Fill(correctp*sin(tt),weight);          // ("fhRecoBinPt","Reconstructed Pt Distribution",100,0,1);
       }
   }// loop over tracks
 
@@ -452,19 +518,23 @@ if (numappjet > 1)
                  alpha = TMath::ACos(crp*ctp*srt*stt+srp*stp*srt*stt+crt*ctt);   
          }
          Double_t correctp = pt[iT]/stt;
-         fhJetPL2->Fill( correctp*cos(alpha));      
+         fhJetPL2->Fill( correctp*cos(alpha),weight);      
          if ( (jethighest->Eta()-eta[iT])*(jethighest->Eta()-eta[iT]) +  
                          (jethighest->Phi()-phi[iT])*(jethighest->Phi()-phi[iT]) < 0.2*0.2 )   
-                 fhJetJT2->Fill( correctp*sin(alpha));
+                 fhJetJT2->Fill( correctp*sin(alpha),weight);
+         fhJetPT2->Fill(correctp*sin(tt),weight);
          if (fNominalEnergy==0.0){
-                 fhFragmFcn2->Fill(  correctp*sin(tt)/parton->Energy()  );
+                 fhFragmFcn2->Fill(  correctp*sin(tt)/(jethighest->Energy()*fScaleFactor),weight  );
          } else
          {
-                  fhFragmFcn2->Fill(  correctp*sin(tt)/fNominalEnergy );
+                  fhFragmFcn2->Fill(  correctp*sin(tt)/fNominalEnergy,weight );
          }
          if (doesJetMeetBinCriteria)
          {
-                 fhRecoBinPt->Fill(correctp*sin(tt));          // ("fhRecoBinPt","Reconstructed Pt Distribution",100,0,1);
+                 fhRecoBinPt->Fill(correctp*sin(tt),weight);          // ("fhRecoBinPt","Reconstructed Pt Distribution",100,0,1);
+                 fhRecoBinPtNoBg->Fill(correctp*sin(tt),weight);          // ("fhRecoBinPt","Reconstructed Pt Distribution",100,0,1);
+                 fhRecoBinFragmFcn->Fill(  correctp*sin(tt)/(jethighest->Energy()*fScaleFactor),weight  ); // This is the jet fragmentation function
+                 fhRecoBinFragmFcnNoBg->Fill(  correctp*sin(tt)/(jethighest->Energy()*fScaleFactor),weight  ); // This is the jet fragmentation function
          }
   }// loop over tracks
   }
@@ -478,18 +548,20 @@ if (numappjet > 1)
   AliEMCALParton* parton;
   AliEMCALJet* jet; 
   jet = jethighest;//fOutput->GetJet(0);
-  fhJetEt->Fill(jet->Energy());
-  fhJetEta->Fill(jet->Eta()  );
-  fhJetPhi->Fill(jet->Phi() );
+  fhJetEt->Fill(jet->Energy()*fScaleFactor,weight);
+  fhJetEtDiff->Fill(jethighest->Energy()*fScaleFactor-et,weight);
+  fhJetInvE->Fill(1.0/(jethighest->Energy()*fScaleFactor),weight);
+  fhJetEta->Fill(jet->Eta(),weight  );
+  fhJetPhi->Fill(jet->Phi(),weight );
   if (nPartons ==0) return;
   parton = fOutput->GetParton(0);
   
-  fhPartonEta->Fill( parton->Eta() );
-  fhPartonPhi->Fill( parton->Phi() );
+  fhPartonEta->Fill( parton->Eta(),weight );
+  fhPartonPhi->Fill( parton->Phi(),weight );
 
   //hJetEtDiff->Fill( jet->Energy() - parton->Energy() );
-  fhEtaDiff->Fill( jet->Eta() - parton->Eta() );
-  fhPhiDiff->Fill( jet->Phi() - parton->Phi() );
+  fhEtaDiff->Fill( jet->Eta() - parton->Eta(),weight );
+  fhPhiDiff->Fill( jet->Phi() - parton->Phi(),weight );
   fhEtaPhiSpread->Fill(jet->Eta()-parton->Eta(),jet->Phi() - parton->Phi());
  /* 
   Float_t *pt,*phi,*eta;
@@ -530,19 +602,20 @@ if (numappjet > 1)
          {
                  alpha = TMath::ACos(crp*ctp*srt*stt+srp*stp*srt*stt+crt*ctt);           }
       Double_t correctp = pt[iT]/stt;  
-      fhPartonPL->Fill( correctp*cos(alpha));
+      fhPartonPL->Fill( correctp*cos(alpha),weight);
       if ( (parton->Eta()-eta[iT])*(parton->Eta()-eta[iT]) +  
       (parton->Phi()-phi[iT])*(parton->Phi()-phi[iT]) < 0.2*0.2 ) 
-             fhPartonJT->Fill( correctp*sin(alpha));
+             fhPartonJT->Fill( correctp*sin(alpha),weight);
       if (fNominalEnergy == 0.0) {
-             fhPartonFragmFcn->Fill(  correctp*sin(tt)/parton->Energy()  );
+             fhPartonFragmFcn->Fill(  correctp*sin(tt)/parton->Energy(),weight  );
+         fhPartonPT->Fill(correctp*sin(tt),weight);
       }else 
       {
-              fhPartonFragmFcn->Fill(correctp*sin(tt)/fNominalEnergy);
+              fhPartonFragmFcn->Fill(correctp*sin(tt)/fNominalEnergy,weight);
       }
       if (doesJetMeetBinCriteria)
       {
-             fhRecoBinPartonPt->Fill(correctp*sin(tt));          // ("fhRecoBinPt","Reconstructed Pt Distribution",100,0,1);
+             fhRecoBinPartonPt->Fill(correctp*sin(tt),weight);          // ("fhRecoBinPt","Reconstructed Pt Distribution",100,0,1);
       }
   }// loop over tracks
 
@@ -575,22 +648,37 @@ if (numappjet > 1)
                  alpha = TMath::ACos(crp*ctp*srt*stt+srp*stp*srt*stt+crt*ctt);   
          }
          Double_t correctp = pt[iT]/stt;
-         fhJetPL->Fill( correctp*cos(alpha));      
+         fhJetPL->Fill( correctp*cos(alpha),weight);      
          if ( (jet->Eta()-eta[iT])*(jet->Eta()-eta[iT]) +  
-                         (jet->Phi()-phi[iT])*(jet->Phi()-phi[iT]) < 0.2*0.2 )   
-                 fhJetJT->Fill( correctp*sin(alpha));
+              (jet->Phi()-phi[iT])*(jet->Phi()-phi[iT]) < 0.2*0.2 )   
+                 fhJetJT->Fill( correctp*sin(alpha),weight);
+         fhJetPT->Fill(correctp*sin(tt),weight);
          if (fNominalEnergy==0.0){
-                 fhFragmFcn->Fill(  correctp*sin(tt)/parton->Energy()  );
+                 fhFragmFcn->Fill(  correctp*sin(tt)/(jet->Energy()*fScaleFactor),weight  ); // This is the jet fragmentation function
          } else
          {
-                  fhFragmFcn->Fill(  correctp*sin(tt)/fNominalEnergy );
+                  fhFragmFcn->Fill(  correctp*sin(tt)/fNominalEnergy,weight );
          }
          if (doesJetMeetBinCriteria)
          {
-                 fhRecoBinPt->Fill(correctp*sin(tt));          // ("fhRecoBinPt","Reconstructed Pt Distribution",100,0,1);
+                 fhRecoBinPt->Fill(correctp*sin(tt),weight);          // ("fhRecoBinPt","Reconstructed Pt Distribution",100,0,1);
+                 fhRecoBinPtNoBg->Fill(correctp*sin(tt),weight);          // ("fhRecoBinPt","Reconstructed Pt Distribution",100,0,1);
+                 fhRecoBinFragmFcn->Fill(  correctp*sin(tt)/(jet->Energy()*fScaleFactor),weight  ); // This is the jet fragmentation function
+                 fhRecoBinFragmFcnNoBg->Fill(  correctp*sin(tt)/(jet->Energy()*fScaleFactor),weight  ); // This is the jet fragmentation function
          }
   }// loop over tracks
   }
+
+if (numappjet>=1 && fhBackHisto != 0 && doesJetMeetBinCriteria)
+  {
+    for (Int_t count=1;count<=100;count++)
+    {
+       fhRecoBinFragmFcnNoBg->Fill( ((Float_t)count)/(jethighest->Energy()*fScaleFactor),-fhBackHisto->GetBinContent(count)*weight);
+       fhRecoBinPtNoBg->AddBinContent(count,-fhBackHisto->GetBinContent(count)*weight);
+    } 
+  }
+
+
 }