Merge pull request cms-sw#47 from amarini/topic_sfth2f

th2f sf. Leftovers

interface/BTagCalibrationStandalone.hpp

@@ -0,0 +1,181 @@
+#ifndef BTagEntry_H
+#define BTagEntry_H
+
+/**
+ *
+ * BTagEntry
+ *
+ * Represents one pt- or discriminator-dependent calibration function.
+ *
+ * measurement_type:    e.g. comb, ttbar, di-mu, boosted, ...
+ * sys_type:            e.g. central, plus, minus, plus_JEC, plus_JER, ...
+ *
+ * Everything is converted into a function, as it is easiest to store it in a
+ * txt or json file.
+ *
+ ************************************************************/
+
+#include <string>
+#include <TF1.h>
+#include <TH1.h>
+
+
+class BTagEntry
+{
+public:
+  enum OperatingPoint {
+    OP_LOOSE=0,
+    OP_MEDIUM=1,
+    OP_TIGHT=2,
+    OP_RESHAPING=3,
+  };
+  enum JetFlavor {
+    FLAV_B=0,
+    FLAV_C=1,
+    FLAV_UDSG=2,
+  };
+  struct Parameters {
+    OperatingPoint operatingPoint;
+    std::string measurementType;
+    std::string sysType;
+    JetFlavor jetFlavor;
+    float etaMin;
+    float etaMax;
+    float ptMin;
+    float ptMax;
+    float discrMin;
+    float discrMax;
+
+    // default constructor
+    Parameters(
+      OperatingPoint op=OP_TIGHT,
+      std::string measurement_type="comb",
+      std::string sys_type="central",
+      JetFlavor jf=FLAV_B,
+      float eta_min=-99999.,
+      float eta_max=99999.,
+      float pt_min=0.,
+      float pt_max=99999.,
+      float discr_min=0.,
+      float discr_max=99999.
+    );
+
+  };
+
+  BTagEntry() {}
+  BTagEntry(const std::string &csvLine);
+  BTagEntry(const std::string &func, Parameters p);
+  BTagEntry(const TF1* func, Parameters p);
+  BTagEntry(const TH1* histo, Parameters p);
+  ~BTagEntry() {}
+  static std::string makeCSVHeader();
+  std::string makeCSVLine() const;
+  static std::string trimStr(std::string str);
+
+  // public, no getters needed
+  std::string formula;
+  Parameters params;
+
+};
+
+#endif  // BTagEntry_H
+
+
+#ifndef BTagCalibration_H
+#define BTagCalibration_H
+
+/**
+ * BTagCalibration
+ *
+ * The 'hierarchy' of stored information is this:
+ * - by tagger (BTagCalibration)
+ *   - by operating point or reshape bin
+ *     - by jet parton flavor
+ *       - by type of measurement
+ *         - by systematic
+ *           - by eta bin
+ *             - as 1D-function dependent of pt or discriminant
+ *
+ ************************************************************/
+
+#include <map>
+#include <vector>
+#include <string>
+#include <istream>
+#include <ostream>
+
+
+class BTagCalibration
+{
+public:
+  BTagCalibration() {}
+  BTagCalibration(const std::string &tagger);
+  BTagCalibration(const std::string &tagger, const std::string &filename);
+  ~BTagCalibration() {}
+
+  std::string tagger() const {return tagger_;}
+
+  void addEntry(const BTagEntry &entry);
+  const std::vector<BTagEntry>& getEntries(const BTagEntry::Parameters &par) const;
+
+  void readCSV(std::istream &s);
+  void readCSV(const std::string &s);
+  void makeCSV(std::ostream &s) const;
+  std::string makeCSV() const;
+
+protected:
+  static std::string token(const BTagEntry::Parameters &par);
+
+  std::string tagger_;
+  std::map<std::string, std::vector<BTagEntry> > data_;
+
+};
+
+#endif  // BTagCalibration_H
+
+
+#ifndef BTagCalibrationReader_H
+#define BTagCalibrationReader_H
+
+/**
+ * BTagCalibrationReader
+ *
+ * Helper class to pull out a specific set of BTagEntry's out of a
+ * BTagCalibration. TF1 functions are set up at initialization time.
+ *
+ ************************************************************/
+
+#include <memory>
+#include <string>
+
+
+
+class BTagCalibrationReader
+{
+public:
+  BTagCalibrationReader() {}
+  BTagCalibrationReader(BTagEntry::OperatingPoint op,
+                        std::string sysType="central");
+
+  void load(const BTagCalibration & c,
+            BTagEntry::JetFlavor jf,
+            std::string measurementType="comb");
+
+  double eval(BTagEntry::JetFlavor jf,
+              float eta,
+              float pt,
+              float discr=0.) const;
+
+  std::pair<float, float> min_max_pt(BTagEntry::JetFlavor jf, 
+                                     float eta, 
+                                     float discr=0.) const;
+
+protected:
+  class BTagCalibrationReaderImpl;
+  std::auto_ptr<BTagCalibrationReaderImpl> pimpl;
+};
+
+
+#endif  // BTagCalibrationReader_H
+
+

interface/Loop.hpp

@@ -100,6 +100,7 @@ class Looper{
         {event_ -> GetWeight()  -> AddPtEtaSF(label,pt1,pt2,eta1,eta2,sf,err); }
         inline void AddSplineSF(string label, double pt, double sf, double err){ event_->GetWeight() -> AddSplineSF(label,pt,sf,err);}
         inline void AddCSVSF(string label, string filename){event_->GetWeight() ->AddCSVSF(label,filename);}
+        inline void AddTh2fSF(string label, string filename){event_->GetWeight() ->AddTh2fSF(label,filename);}
 
         // -- PU Reweight
         inline void AddTarget( TH1*h, int runMin=-1, int runMax =-1,double lumi=-1){ event_ -> GetWeight() -> AddTarget(h,runMin,runMax,lumi);}

interface/RooCMSShape.hpp

@@ -0,0 +1,61 @@
+/*****************************************************************************
+ * Project: CMS detector at the CERN
+ *
+ * Package: PhysicsTools/TagAndProbe/RooCMSShape
+ *
+ *
+ * Authors:
+ *   Nadia Adam, Princeton - neadam@princeton.edu
+ *   Adam Hunt, Princeton  - ahunt@princeton.edu
+ *   Kalanand Mishra, Fermilab - kalanand@fnal.gov
+ *
+ * Description:
+ *   Defines a probability density function which has exponential decay 
+ *   distribution at high mass beyond the pole position (say, Z peak)  
+ *   but turns over (i.e., error function) at low mass due to threshold 
+ *   effect. We use this to model the background shape in Z->ll invariant 
+ *   mass.
+ * History:
+ *   
+ *
+ * Copyright (C) 2008 FNAL 
+ *****************************************************************************/
+
+#ifndef ROO_CMS_SHAPE
+#define ROO_CMS_SHAPE
+
+#include "RooAbsPdf.h"
+#include "RooRealProxy.h"
+#include "RooAbsReal.h"
+#include "TMath.h"
+#include "RooMath.h"
+
+class RooCMSShape : public RooAbsPdf {
+public:
+  RooCMSShape() {}
+  RooCMSShape(const char *name, const char *title,
+	      RooAbsReal& _x,
+	      RooAbsReal& _alpha,
+	      RooAbsReal& _beta,
+	      RooAbsReal& _gamma,
+	      RooAbsReal& _peak);
+
+  RooCMSShape(const RooCMSShape& other, const char* name);
+  inline virtual TObject* clone(const char* newname) const { return new RooCMSShape(*this,newname); }
+  inline ~RooCMSShape() {}
+  Double_t evaluate() const ;
+
+
+//  ClassDef(RooCMSShape,1);
+
+protected:
+
+  RooRealProxy x ;
+  RooRealProxy alpha ;
+  RooRealProxy beta ;
+  RooRealProxy gamma ;
+  RooRealProxy peak ;
+
+};
+
+#endif

interface/RooVoigtianShape.hpp

@@ -0,0 +1,53 @@
+/*****************************************************************************
+ * Project: RooFit                                                           *
+ * Package: RooFitModels                                                     *
+ *    File: $Id: RooCBShape.h,v 1.11 2007/07/12 20:30:49 wouter Exp $
+ * Authors:                                                                  *
+ *   WV, Wouter Verkerke, UC Santa Barbara, verkerke@slac.stanford.edu       *
+ *   DK, David Kirkby,    UC Irvine,         dkirkby@uci.edu                 *
+ *                                                                           *
+ * Copyright (c) 2000-2005, Regents of the University of California          *
+ *                          and Stanford University. All rights reserved.    *
+ *                                                                           *
+ * Redistribution and use in source and binary forms,                        *
+ * with or without modification, are permitted according to the terms        *
+ * listed in LICENSE (http://roofit.sourceforge.net/license.txt)             *
+ *****************************************************************************/
+#ifndef ROO_VOIGTIAN_SHAPE
+#define ROO_VOIGTIAN_SHAPE
+
+#include "RooAbsPdf.h"
+#include "RooRealProxy.h"
+
+class RooRealVar;
+
+class RooVoigtianShape : public RooAbsPdf {
+public:
+  RooVoigtianShape() {}
+  RooVoigtianShape(const char *name, const char *title, RooAbsReal& _m,
+	     RooAbsReal& _m0, RooAbsReal& _sigma,
+		   RooAbsReal& _alpha, RooAbsReal& _n,RooAbsReal& _width,Bool_t doFast);
+
+  RooVoigtianShape(const RooVoigtianShape& other, const char* name = 0);
+  virtual TObject* clone(const char* newname) const { return new RooVoigtianShape(*this,newname); }
+
+  inline virtual ~RooVoigtianShape() { }
+
+protected:
+
+
+  RooRealProxy m;
+  RooRealProxy m0;
+  RooRealProxy sigma;
+  RooRealProxy alpha;
+  RooRealProxy n;
+  RooRealProxy width;
+  Double_t evaluate() const;
+  Double_t voigtian(Double_t iX) const;
+private:
+  Bool_t   _doFast;
+  Double_t _invRootPi;
+//  ClassDef(RooVoigtianShape,1) // Crystal Ball lineshape PDF
+};
+
+#endif

interface/SF.hpp

@@ -56,6 +56,16 @@ class SF_PtEta : virtual public SF
 
 const bool operator<( const SF_PtEta::range&r1 , const SF_PtEta::range &r2);
 
+
+class SF_TH2F : virtual public SF_PtEta
+{
+    public:
+        SF_TH2F() : SF_PtEta(),SF(){}
+        SF_TH2F(string filename) : SF_PtEta(),SF(){ init(filename); }
+        void init(string filename,string histname="EGamma_SF2D");
+        const string name() const {return "SF_TH2F";}
+};
+
 #include "TSpline.h"
 
 class SF_PtSpline : virtual public SF

interface/Weight.hpp

@@ -47,6 +47,7 @@ class Weight{
 
     void AddSF( string label, double sf, double err);
     void AddPtEtaSF( string label, double pt1,double pt2 ,double eta1 ,double eta2,double sf, double err);
+    void AddTh2fSF(string label, string filename);
     void AddSplineSF(string label, double pt, double sf, double err);
     void AddCSVSF(string label, string filename);
 

python/Loop.py

@@ -181,6 +181,9 @@ def applyJson(obj,fname):
 	if key['type'] == 'pteta':
 		if opts.verbose: print label,key['type'],  key['pt1'] ,key['pt2'],key['eta1'],key['eta2'],key['sf'],key['err'] 
 		loop.AddPtEtaSF(label, key['pt1'] ,key['pt2'],key['eta1'],key['eta2'],key['sf'],key['err'])
+	if key['type'] == 'th2f':
+		if opts.verbose: print label,key['type'],  key['filename']
+		loop.AddTh2fSF(label, key['filename'])
 	if key['type'] == 'base':
 		if opts.verbose: print label,key['type'], key['sf'],key['err'] 
 		loop.AddSF(label, key['sf'], key['err'])

python/ParseDat.py

@@ -287,6 +287,11 @@ def ReadSFDB(file):
 			err  = float ( l.split(' ')[4] )
 			R['pt']=pt
 
+		elif type == 'th2f':
+			R['filename'] = l.split(' ' )[2]
+			sf=0.0 ## ignored
+			err=0.0 ## ignored
+
 		elif type == 'base':
 			sf  = float ( l.split(' ') [2] )
 			err = float ( l.split(' ') [3] )

script/doTagAndProbe.py

@@ -148,7 +148,7 @@ def HistName( type, mc ,pt , Pass=True):
 		c2=ROOT.TCanvas("c2","c")
 		data = ROOT.RooDataHist("h","h",ROOT.RooArgList(x),hists[name])
 		#bkg = ROOT.RooCMSShape("cms-shape","cms-shape",x,Balpha,Bbeta,Bgamma,Bpeak)
-		bkg = ROOT.RooBerneinFast(2)("bern","bern",x,ROOT.RooArgList(Bbern0,Bbern1,Bbern2))
+		bkg = ROOT.RooBernsteinFast(2)("bern","bern",x,ROOT.RooArgList(Bbern0,Bbern1,Bbern2))
 		# RooVoigtianShape::RooVoigtianShape(const char* name, const char* title, RooAbsReal& _m, RooAbsReal& _m0, RooAbsReal& _sigma, RooAbsReal& _alpha, RooAbsReal& _n, RooAbsReal& _width, bool doFast) =>
 		sig = ROOT.RooVoigtianShape("voigtian","voigtian", x,Zmass, Ssigma, Salpha,Sn,Zwidth,True)
 		f.setVal(.9)