ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/cbrown/AnalysisFramework/Plotting/Modules/Systematics.C
(Generate patch)

Comparing UserCode/cbrown/AnalysisFramework/Plotting/Modules/Systematics.C (file contents):
Revision 1.46 by buchmann, Wed Nov 9 10:55:25 2011 UTC vs.
Revision 1.51 by buchmann, Thu Nov 10 15:16:32 2011 UTC

# Line 149 | Line 149 | float allcontributionsplot(TTree* events
149   TH1F* plotEff(TTree* events, TCut kbase, TString informalname, int flipped) {
150          iplot++;
151          int count=iplot;
152 +        iplot++;
153 +        int count2=iplot;
154          // Define new histogram
155          char hname[30]; sprintf(hname,"hJzbEff%d",count);
156 <        TH1F* hJzbEff = new TH1F(hname,"JZB selection efficiency ; JZB (GeV/c); Efficiency",
157 <                                                         nBins,jzbMin,jzbMax);
156 >        char hname2[30]; sprintf(hname2,"hJzbEff%d",count2);
157 >        TH1F* hJzbEff = new TH1F(hname,"JZB selection efficiency ; JZB (GeV/c); Efficiency",nBins,jzbMin,jzbMax);
158 >        TH1F* hJzbEff2= new TH1F(hname2,"JZB selection efficiency ; JZB (GeV/c); Efficiency",1,-14000,14000);
159          Float_t step = (jzbMax-jzbMin)/static_cast<Float_t>(nBins);
160 <        
161 <        if(flipped==0) events->Draw(mcjzbexpression.c_str(),"genJZB>-400"&&kbase,"goff");
162 <        else events->Draw(("-"+mcjzbexpression).c_str(),"genJZB>-14000"&&kbase,"goff");
163 <        Float_t maxEff = events->GetSelectedRows();
160 >
161 >        if(flipped==0) events->Draw((mcjzbexpression+">>"+(string)hname2).c_str(),("genJZB>-400"&&kbase),"goff");
162 >        else events->Draw(("(-"+mcjzbexpression+")>>"+(string)hname2).c_str(),("genJZB>-400"&&kbase),"goff");
163 >        Float_t maxEff = hJzbEff2->Integral();
164          if(verbose>0) dout << hname << " (" << informalname <<") " << maxEff <<  std::endl;
165          
166          if(verbose>0) dout <<  "JZB max = " << jzbMax << std::endl;
# Line 165 | Line 168 | TH1F* plotEff(TTree* events, TCut kbase,
168          char cut[256];
169          for ( Int_t iBin = 0; iBin<nBins; ++iBin ) {
170                  sprintf(cut,"genJZB>%3f",jzbMin+iBin*step);
171 <                events->Draw(mcjzbexpression.c_str(),TCut(cut)&&kbase,"goff");
172 <                Float_t eff = static_cast<Float_t>(events->GetSelectedRows())/maxEff;
173 <                //     dout << "COUCOU " << __LINE__ << std::endl;
171 >                if(flipped==0) events->Draw((mcjzbexpression+">>"+(string)hname2).c_str(),(TCut(cut)&&kbase),"goff");
172 >                if(flipped>0) events->Draw(("(-"+mcjzbexpression+")>>"+(string)hname2).c_str(),(TCut(cut)&&kbase),"goff");
173 >                Float_t eff = static_cast<Float_t>(hJzbEff2->Integral())/maxEff;
174                  hJzbEff->SetBinContent(iBin+1,eff);
175                  hJzbEff->SetBinError(iBin+1,TMath::Sqrt(eff*(1-eff)/maxEff));
176          }
177 +        delete hJzbEff2;
178          return hJzbEff;
175 –        
176 –        
179   }
180  
181  
# Line 252 | Line 254 | float pileup(TTree *events, bool require
254          
255          // Pimped-up function
256          TF1* funcUp = (TF1*)func->Clone();
257 <        funcUp->SetParameter( 0., func->GetParameter(0)/1.1); // 10% systematic error (up in sigma => 0.1 in erfc)
257 >        funcUp->SetParameter( 0, func->GetParameter(0)/1.1); // 10% systematic error (up in sigma => 0.1 in erfc)
258          if(!automatized) dout << "  PU: " << funcUp->Eval(jzbSel) << " " <<  func->Eval(jzbSel)
259          << "(" << (funcUp->Eval(jzbSel)-func->Eval(jzbSel))/func->Eval(jzbSel)*100. << "%)" << std::endl;
260          
# Line 294 | Line 296 | void PeakError(TTree *events,float &resu
296          }
297          if(TMath::Abs(rescent-resup)>TMath::Abs(rescent-resdown)) result=(TMath::Abs(rescent-resup)/(float)rescent);
298          else result=(TMath::Abs(rescent-resdown)/(float)rescent);
299 +        cout << "    " << result << endl;
300   }
301  
302  
# Line 424 | Line 427 | Value MCefficiency(TTree *events,float &
427          string emnewNegSide = "((id1!=id2)&&(" + snegSide + "))*" + svar; // only used for off peak analysis
428  
429          TH1F *effh= new TH1F("effh","effh",1,-14000,14000);
430 <        if(k>=0)events->Draw((mcjzbexpression+">>effh").c_str(), newPosSide.c_str(),"goff");
431 <        else events->Draw((mcjzbexpression+">>effh").c_str(), (sposSide+"&&(id1==id2)").c_str(),"goff");//the OSSF condition is added for the offpeak analysis, in onpeak case it's there already but doesn't change anything.
430 >        if(k>=0)events->Draw((mcjzbexpression+">>effh").c_str(), TCut(newPosSide.c_str())*PlottingSetup::Weight,"goff");
431 >        else events->Draw((mcjzbexpression+">>effh").c_str(), TCut((sposSide+"&&(id1==id2)").c_str())*PlottingSetup::Weight,"goff");//the OSSF condition is added for the offpeak analysis, in onpeak case it's there already but doesn't change anything.
432          Float_t sel = effh->Integral();
433          Float_t nsel=0;
434          
# Line 434 | Line 437 | Value MCefficiency(TTree *events,float &
437          if(ConsiderSignalContaminationForLimits) {
438            flag_this_change(__FUNCTION__,__LINE__,true);//PlottingSetup::RestrictToMassPeak
439            if(PlottingSetup::RestrictToMassPeak) {
440 <            events->Draw((mcjzbexpression+">>effh").c_str(), newNegSide.c_str(),"goff");
440 >            events->Draw((mcjzbexpression+">>effh").c_str(), (TCut(newNegSide.c_str()))*PlottingSetup::Weight,"goff");
441              nsel += effh->Integral();
442            } else {
443 <            events->Draw((mcjzbexpression+">>effh").c_str(), newNegSide.c_str(),"goff");
443 >            events->Draw((mcjzbexpression+">>effh").c_str(), (TCut(newNegSide.c_str()))*PlottingSetup::Weight,"goff");
444              nsel += effh->Integral();
445 <            events->Draw((mcjzbexpression+">>effh").c_str(), emnewPosSide.c_str(),"goff");
445 >            events->Draw((mcjzbexpression+">>effh").c_str(), (TCut(emnewPosSide.c_str()))*PlottingSetup::Weight,"goff");
446              nsel += effh->Integral();
447 <            events->Draw((mcjzbexpression+">>effh").c_str(), emnewNegSide.c_str(),"goff");
447 >            events->Draw((mcjzbexpression+">>effh").c_str(), (TCut(emnewNegSide.c_str()))*PlottingSetup::Weight,"goff");
448              nsel -= effh->Integral();
449            }
450          }
# Line 506 | Line 509 | void JZBefficiency(TTree *events, string
509  
510   //________________________________________________________________________
511   // Effect of energy scale on efficiency
512 < void JZBjetScale(TTree *events, float &jesdown, float &jesup, string informalname, int flipped, bool requireZ,string addcut="",float syst=0.1, Float_t jzbSelection=-1, TString plotName = "" ) {
512 > void JZBjetScale(TTree *events, float &jesdown, float &jesup, string informalname, int flipped, bool requireZ,string addcut="",Float_t jzbSelection=-1, TString plotName = "" ) {
513          TCut kbase(genMassCut&&"genZPt>0");
514          if(addcut!="") kbase=kbase&&addcut.c_str();//this is mostly for SUSY scans (adding requirements on masses)
515          flag_this_change(__FUNCTION__,__LINE__,true);//PlottingSetup::RestrictToMassPeak
# Line 515 | Line 518 | void JZBjetScale(TTree *events, float &j
518          TCut ksel(cutmass&&cutOSSF);
519          TCut nJets("pfJetGoodNum>2");
520          stringstream down,up;
521 <        down << "pfJetGoodNum"<<30*(1-syst)<<">=3";
522 <        up << "pfJetGoodNum"<<30*(1+syst)<<">=3";
521 >        down << "pfJetGoodNumn1sigma>=3";
522 >        up << "pfJetGoodNump1sigma>=3";
523          
524          TCut nJetsP(up.str().c_str());
525          TCut nJetsM(down.str().c_str());
# Line 524 | Line 527 | void JZBjetScale(TTree *events, float &j
527          if ( !(plotName.Length()>1) ) plotName = informalname;
528          
529          nBins = 1; jzbMin = jzbSel*0.95; jzbMax = jzbSel*1.05;
530 <        TH1F* hist = plotEff(events,(kbase&&ksel&&nJets),informalname,flipped);
528 <        
530 >        TH1F* hist  = plotEff(events,(kbase&&ksel&&nJets),informalname,flipped);
531          TH1F* histp = plotEff(events,(kbase&&ksel&&nJetsP),informalname,flipped);
530 –        
532          TH1F* histm = plotEff(events,(kbase&&ksel&&nJetsM),informalname,flipped);
533          
534          // Dump some information
# Line 637 | Line 638 | int get_npdfs(TTree *events) {
638    
639  
640   void do_systematics_for_one_file(TTree *events,int Neventsinfile,string informalname, vector<vector<float> > &results,int flipped, string mcjzb,string datajzb,float peakerror,bool requireZ=false, string addcut="", bool ismSUGRA=false) {
641 <  float JetEnergyScaleUncert=0.1;
641 <  float JZBScaleUncert=0.1;
641 >  float JZBScaleUncert=0.05;
642    mcjzbexpression=mcjzb;
643 <  float triggereff=5.0/100;// in range [0,1]
643 >  float triggereff=2.0/100;// in range [0,1]
644    dout << "Trigger efficiency not implemented in this script  yet, still using external one" << endl;
645    float leptonseleff=2.0/100;// in range [0,1]
646    leptonseleff=TMath::Sqrt(leptonseleff*leptonseleff+leptonseleff*leptonseleff); // because the 2% is per lepton
# Line 658 | Line 658 | void do_systematics_for_one_file(TTree *
658    if(PlottingSetup::computeJZBefficiency) JZBefficiency(events,informalname,jzbeff,jzbefferr,flipped,requireZ,addcut);
659    if(!automatized) dout << "JZB efficiency: " << jzbeff << "+/-" << jzbefferr << endl;
660    
661 <  if(!automatized) dout << "Error from Peak position:" << endl;
661 >  if(!automatized) dout << "Error from Peak position:";
662    float sysfrompeak=0;
663    PeakError(events,sysfrompeak,mcjzb,peakerror,flipped,addcut);
664      
665 <  if(!automatized) dout << "Jet energy scale: " << std::endl;
665 >  if(!automatized) dout << "Jet energy scale (JES): " << std::endl;
666    float jesup,jesdown;
667 <  JZBjetScale(events,jesdown,jesup,informalname,flipped,requireZ,addcut,JetEnergyScaleUncert);
667 >  JZBjetScale(events,jesdown,jesup,informalname,flipped,requireZ,addcut);
668    
669    if(!automatized) dout << "JZB scale: " << std::endl;
670    float scaleup,scaledown,scalesyst;
# Line 674 | Line 674 | void do_systematics_for_one_file(TTree *
674    float resp,resperr;
675    if(PlottingSetup::computeJZBresponse) {
676          if(!automatized) dout << "JZB response: " << std::endl;
677 <        JZBresponse(events,requireZ,resp,resperr,flipped,addcut);
677 >        if(!ismSUGRA) JZBresponse(events,requireZ,resp,resperr,flipped,addcut);
678    }
679  
680    if(!automatized) dout << "Pileup: " << std::endl;
681 <  float resolution;
682 <  resolution=pileup(events,requireZ,informalname,flipped,addcut);
681 > //  float resolution;
682 >  //resolution=pileup(events,requireZ,informalname,flipped,addcut);
683  
684    float PDFuncert=0;
685    if(!automatized) dout << "Assessing PDF uncertainty: " << std::endl;
# Line 694 | Line 694 | void do_systematics_for_one_file(TTree *
694    dout << "Lepton Sel Eff: " << leptonseleff << endl; // in range [0,1]
695    dout << "Jet energy scale: " << jesup << " " << jesdown << endl; // in range [0,1]
696    dout << "JZB Scale Uncert: " << scaledown << " " << scaleup << endl; // in range [0,1]
697 <  dout << "Resolution : " << resolution << endl; // in range [0,1]
697 > //  dout << "Resolution : " << resolution << endl; // in range [0,1]
698    dout << "From peak : " << sysfrompeak << endl; // in range [0,1]
699    if(ismSUGRA) dout << "PDF uncertainty  : " << PDFuncert << endl; // in range [0,1]
700    if(PlottingSetup::computeJZBefficiency) dout << "JZB efficiency: " << jzbeff << "+/-" << jzbefferr << " (not yet included below) " << endl; // in range [0,1]
# Line 705 | Line 705 | void do_systematics_for_one_file(TTree *
705    toterr+=(leptonseleff)*(leptonseleff);
706    if(fabs(jesup)>fabs(jesdown)) toterr+=(jesup*jesup); else toterr+=(jesdown*jesdown);
707    if(fabs(scaleup)>fabs(scaledown)) toterr+=(scaleup*scaleup); else toterr+=(scaledown*scaledown);
708 <  toterr+=(resolution*resolution);
708 > //  toterr+=(resolution*resolution);
709    toterr+=(sysfrompeak*sysfrompeak);
710    if(ismSUGRA) toterr+=(PDFuncert*PDFuncert);
711    dout << "TOTAL SYSTEMATICS: " << TMath::Sqrt(toterr) << " --> " << TMath::Sqrt(toterr)*mceff << endl;
# Line 725 | Line 725 | void do_systematics_for_one_file(TTree *
725    res.push_back(TMath::Sqrt((mcefferr)*(mcefferr)+(toterr*toterr)));
726    if(fabs(jesup)>fabs(jesdown)) res.push_back(fabs(jesup)); else res.push_back(fabs(jesdown));
727    if(fabs(scaleup)>fabs(scaledown)) res.push_back(fabs(scaleup)); else res.push_back(fabs(scaledown));
728 <  res.push_back(fabs(resolution));
728 > //  res.push_back(fabs(resolution));
729 >  res.push_back(0.0);
730    res.push_back(mceff_nosigcont.getValue());
731    res.push_back(mceff_nosigcont.getError());
732    if(ismSUGRA) res.push_back(PDFuncert);

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines