Resync with Wolfgang's fork

DrawHits/drawHits.py

@@ -254,6 +254,8 @@ def fillHistoByDef(tree,hDef,extraCuts):
         if hDef.vetoMType(mType):
             continue
         is1D = hDef.getParameter('yNbins',mType)==None
+        is2D = hDef.getParameter('yNbins',mType)!=None and hDef.getParameter('zNbins',mType)==None
+        is3D = hDef.getParameter('yNbins',mType)!=None and hDef.getParameter('zNbins',mType)!=None
         isProfile = hDef.getParameter('profile',mType)!=None and hDef.getParameter('profile',mType)
         effCuts = hDef.getParameter('effCuts',mType)
         variable = hDef.getParameter('variable',mType)
@@ -278,6 +280,7 @@ def fillHistoByDef(tree,hDef,extraCuts):
                 tree.Project(hName+"_2",variable, \
                             cutString(extraCuts,hDef.getParameter('baseCuts',mType),"moduleType=="+str(mType),effCuts))
                 histos[mType][3] = ROOT.TEfficiency(histos[mType][1],histos[mType][0])
+                histos[mType][3].SetMarkerStyle(20)
             else:
                 # always keep final histogram in 4th position
                 histos[mType][3] = histos[mType][0]
@@ -290,7 +293,7 @@ def fillHistoByDef(tree,hDef,extraCuts):
                           cutString(extraCuts,hDef.getParameter('baseCuts',mType),"moduleType=="+str(mType)))
             # always keep final histogram in 4th position
             histos[mType][3] = histos[mType][0]
-        else:
+        elif is2D:
             nby = hDef.getParameter('yNbins',mType)
             ymin = hDef.getParameter('yMin',mType)
             ymax = hDef.getParameter('yMax',mType)
@@ -310,13 +313,27 @@ def fillHistoByDef(tree,hDef,extraCuts):
             else:
                 # always keep final histogram in 4th position
                 histos[mType][3] = histos[mType][0]
+        elif is3D:
+            nby = hDef.getParameter('yNbins',mType)
+            ymin = hDef.getParameter('yMin',mType)
+            ymax = hDef.getParameter('yMax',mType)
+            nbz = hDef.getParameter('zNbins',mType)
+            zmin = hDef.getParameter('zMin',mType)
+            zmax = hDef.getParameter('zMax',mType)
+            histos[mType] = [ ROOT.TH3F(hName+"_1",hName+"_1",nbx,xmin,xmax,nby,ymin,ymax,nbz,zmin,zmax), \
+                                  None, None, None ]
+            tree.Project(hName+"_1",variable, \
+                          cutString(extraCuts,hDef.getParameter('baseCuts',mType),"moduleType=="+str(mType)))
+            assert effCuts==None
+            # always keep final histogram in 4th position
+            histos[mType][3] = histos[mType][0]
         #print("Ending for ",hDef.name,hName,hTitle)
 
     savedDir.cd()
     return histos
 
 
-def drawHistoByDef(histos,hDef,logY=False,same=False):
+def drawHistoByDef(histos,hDef,logY=False,logZ=False,same=False):
     result = { 'cnv' : None, 'histos' : histos, 'pave' : None }
 
     savedDir = ROOT.gDirectory
@@ -344,6 +361,8 @@ def drawHistoByDef(histos,hDef,logY=False,same=False):
         if hDef.vetoMType(mType):
             continue
         is1D = hDef.getParameter('yNbins',mType)==None
+        is2D = hDef.getParameter('yNbins',mType)!=None and hDef.getParameter('zNbins',mType)==None
+        is3D = hDef.getParameter('yNbins',mType)!=None and hDef.getParameter('zNbins',mType)!=None
         isProfile = hDef.getParameter('profile',mType)!=None and hDef.getParameter('profile',mType)
         effCuts = hDef.getParameter('effCuts',mType)
         variable = hDef.getParameter('variable',mType)
@@ -362,6 +381,7 @@ def drawHistoByDef(histos,hDef,logY=False,same=False):
         xmax = hDef.getParameter('xMax',mType)
 
         ytitle = hDef.getParameter('yTitle',mType) if hDef.getParameter('yTitle',mType) else ""
+        ztitle = hDef.getParameter('zTitle',mType) if hDef.getParameter('zTitle',mType) else ""
         if is1D and ( not isProfile ):
             ymin = hDef.getParameter('yMin',mType) if hDef.getParameter('yMin',mType)!=None else 0.
             ymax = hDef.getParameter('yMax',mType) if hDef.getParameter('yMax',mType)!=None else 1.05
@@ -372,13 +392,16 @@ def drawHistoByDef(histos,hDef,logY=False,same=False):
                 histos[mType][2].GetXaxis().SetTitle(xtitle)
                 histos[mType][2].GetYaxis().SetTitle(ytitle)
                 histos[mType][3] = ROOT.TEfficiency(histos[mType][1],histos[mType][0])
-                histos[mType][3].SetMarkerSize(0.3)
+                histos[mType][3].SetMarkerStyle(20)
                 histos[mType][3].Draw("same Z")
             else:
                 histos[mType][0].SetTitle(hTitle)
                 histos[mType][0].GetXaxis().SetTitle(xtitle)
                 histos[mType][0].GetYaxis().SetTitle(ytitle)
                 histos[mType][0].Draw("same" if same else "")
+            fitFunc = hDef.getParameter('fit')
+            if fitFunc!=None:
+                histos[mType][0].Fit(fitFunc,"Q","same")
         elif isProfile:
             histos[mType][0].SetTitle(hTitle)
             histos[mType][0].GetXaxis().SetTitle(xtitle)
@@ -387,7 +410,7 @@ def drawHistoByDef(histos,hDef,logY=False,same=False):
             #histos[mType][0].SetFillColor(ROOT.TColor.GetColorBright(ROOT.kGray))
             #histos[mType][0].SetFillColor(ROOT.kGray)
             histos[mType][0].Draw("same" if same else "")
-        else:
+        elif is2D:
             assert not same
             zmin = hDef.getParameter('zMin',mType) if hDef.getParameter('zMin',mType)!=None else 0.
             zmax = hDef.getParameter('zMax',mType) if hDef.getParameter('zMax',mType)!=None else 1.05
@@ -403,8 +426,20 @@ def drawHistoByDef(histos,hDef,logY=False,same=False):
                 histos[mType][0].GetXaxis().SetTitle(xtitle)
                 histos[mType][0].GetYaxis().SetTitle(ytitle)
                 histos[mType][0].Draw("ZCOL")
+        elif is3D:
+            assert not same
+            zmin = hDef.getParameter('zMin',mType) if hDef.getParameter('zMin',mType)!=None else 0.
+            zmax = hDef.getParameter('zMax',mType) if hDef.getParameter('zMax',mType)!=None else 1.05
+            assert effCuts==None
+            histos[mType][0].SetTitle(hTitle)
+            histos[mType][0].GetXaxis().SetTitle(xtitle)
+            histos[mType][0].GetYaxis().SetTitle(ytitle)
+            histos[mType][0].GetZaxis().SetTitle(ztitle)
+            histos[mType][0].Draw()
         if logY or hDef.getParameter('logY',mType):
             ROOT.gPad.SetLogy(1)
+        if is2D and ( logZ or hDef.getParameter('logZ',mType) ):
+            ROOT.gPad.SetLogz(1)
         ROOT.gPad.Update()
 
         #cnv.cd()
@@ -484,10 +519,13 @@ def addHistogram(varString,cuts,effCuts=None,name='userHist'):
                         type=str, default='*')
 parser.add_argument('--vetoedHistograms', help='comma-separated names of histogram definitions not to be used',
                         type=str, default='')
-parser.add_argument('--logY', help='use log scale', action='store_true', default=False)
+parser.add_argument('--logY', help='use log scale for y axis', action='store_true', default=False)
+parser.add_argument('--logZ', help='use log scale for z axis', action='store_true', default=False)
 parser.add_argument('--printTree', '-p', help='print TTree contents', action='store_true', default=False)
 parser.add_argument('--listHistograms', '-l', help='list predefined and selected histograms', \
                         action='store_true', default=False)
+parser.add_argument('--zone', '-z', help='restrict to zone in OT (barrel, tilted, endcap)', type=str, \
+                        choices=['barrel','tilted','endcap'], default=None)
 parser.add_argument('file', help='input file', type=str, nargs='+', default=None)
 args = parser.parse_args()
 outputFormats = [ ]
@@ -498,6 +536,19 @@ def addHistogram(varString,cuts,effCuts=None,name='userHist'):
 fitResiduals = args.fitResiduals.split(",") if args.fitResiduals else [ ]
 selectedHistoNames = args.selectedHistograms.split(",")
 vetoedHistoNames = args.vetoedHistograms.split(",")
+#
+# add cut for zone definition?
+#
+if args.zone=="barrel":
+    zoneCuts = "detNormal.Rho()>0.99"
+elif args.zone=="endcap":
+    zoneCuts = "detNormal.Rho()<0.01"
+elif args.zone=="tilted":
+    zoneCuts = "detNormal.Rho()>0.05&&detNormal.Rho()<0.095"
+else:
+    zoneCuts = ""
+args.cuts = cutString(args.cuts,zoneCuts)
+
 #
 # load histogram definitions
 #
@@ -591,14 +642,17 @@ def addHistogram(varString,cuts,effCuts=None,name='userHist'):
     for hName in allHDefs.byCanvas[cName]:
         print("Processing histogram",hName,"in canvas",cName)
         cHistos[hName] = fillHistoByDef(simHitTree,allHDefs.byCanvas[cName][hName],extraCuts)
-        allObjects.append(drawHistoByDef(cHistos[hName],allHDefs.byCanvas[cName][hName],logY=args.logY,same=same))
+        allObjects.append(drawHistoByDef(cHistos[hName],allHDefs.byCanvas[cName][hName], \
+                                             logY=args.logY,logZ=args.logZ,same=same))
         same = True
     if args.output!=None:
         c = allObjects[-1]['cnv']
         #for c in [ x['cnv'] for x in allObjects ]:
         basename = os.path.join(args.output,c.GetName())
         if args.sampleName!=None:
             basename += "_" + args.sampleName
+        if args.zone!=None:
+            basename += "_" + args.zone
         print(basename)
         for fmt in outputFormats:
             c.SaveAs(basename+"."+fmt)

HitAnalyzer/test/histogramDefinition.py → DrawHits/drawHitsConfiguration.py

@@ -1,20 +1,30 @@
 #
 # Definitions for histograms from a configuration file
 #
+import yaml
 from fnmatch import fnmatch
 
 class HistogramDefinition:
     ''' A single histogram definition. 
     '''
-    reqGenFields = [ 'canvasName', 'histogramName', 'histogramTitle', 'variable', 'baseCuts' ]
+    # fields that have to be present in the general section
+    reqGenFields = [ 'canvasName' ]
+    # fields that have to be present in a histogram section
     reqHistFields = [ ]
-    requiredFields = reqGenFields + reqHistFields
-    optGenFields =  [ 'effCuts', 'logY' ]
-    optHistFields = [ 'xNbins', 'xMin', 'xMax', 'xTitle', 'yTitle', 'yNbins', 'yMin', 'yMax', \
-                          'zMin', 'zMax', 'display', 'profile' ]
-    optionalFields = optGenFields + optHistFields
-    allFields = requiredFields + optionalFields
-    allHistFields = reqHistFields + optHistFields
+    requiredFields = list(set(reqGenFields + reqHistFields))
+    # optional fields in the general section
+    optGenFields =  [ 'variable', 'baseCuts', 'effCuts', 'logY', 'logZ' ]
+    # optional fields in the histogram section
+    optHistFields = [ 'variable', 'histogramName', 'histogramTitle', 'fit', \
+                          'xNbins', 'xMin', 'xMax', 'xTitle', 'yTitle', 'zTitle', \
+                          'yNbins', 'yMin', 'yMax', \
+                          'zNbins', 'zMin', 'zMax', 'display', 'profile' ]
+    optionalFields = list(set(optGenFields + optHistFields))
+    allFields = list(set(requiredFields + optionalFields))
+    allHistFields = list(set(reqHistFields + optHistFields))
+    # fields that cannot be present for a single module type
+    vetoMtypeFields = [ 'variable', 'baseCuts', 'effCuts', 'profile' ]
+    vetoMtypeFields = [ 'profile' ]
 
     def __init__(self,name,inputDict):
         ''' Define histogram and drawing parameters from dictionary.
@@ -98,6 +108,9 @@ def getParameter(self,name,mType=None):
         if ( mTName in self.parameters ) and ( name in self.parameters[mTName] ):
           result = self.parameters[mTName][name]
           if result!=None:
+            # check for parameters that can only be general
+            if name in HistogramDefinition.vetoMtypeFields:
+              raise Exception('Parameter '+name+' cannot be present in the section for an individual module type')
             return self.parameters[mTName][name]
         #
         # not found: use general parameter
@@ -106,6 +119,11 @@ def getParameter(self,name,mType=None):
             return self.parameters[name]
         return None
 
+    def __call__(self,name,mType=None):
+        ''' Make access to parameters easier - use function call
+        '''
+        return self.getParameter(name,mType)
+
     def vetoMType(self,mType):
         ''' Check for an mType entry with display = False
         '''
@@ -151,29 +169,124 @@ def __getitem__(self,name):
             return self.allDefinitions[name]
         return None
 
-def loadHistogramDefinitions(configName,selectedNames=[],vetoedNames=[]):
-    ''' Load histogram definitions from a configuration file. The configuration file
+class VarMaskCombinations:
+    ''' Combinations of a variable and a mask for use with RDF.Histo*.
+    '''
+    def __init__(self):
+        #
+        # variable+mask name indexed by variable+selection string
+        #
+        self.varMaskDefinitions = { }
+
+    def __call__(self,rdf,varName,selection):
+        ''' Return the name of a variable+mask combination. Define if necessary.
+            Arguments:
+              rdf ........ RDataFrame to be used for definition
+              varName .... name of the variable
+              selection .. selection string to be used as mask
+        '''
+        #
+        # use variable name and normalized selection string (remove spaces)
+        #
+        selNorm = selection.replace(" ","")
+        varMask = "(" + varName + ")["+selNorm+"]"
+        if varMask in self.varMaskDefinitions:
+            # combination exists (assume that it's also defined in the RDataFrame)
+            varMaskName = self.varMaskDefinitions[varMask]
+            assert varMaskName in rdf.GetDefinedColumnNames()
+        else:
+            # create name for combination and define it in the RDataFrame
+            n = len(self.varMaskDefinitions)
+            varMaskName = "varMask{:03d}".format(n)
+            self.varMaskDefinitions[varMask] = varMaskName
+            rdf = rdf.Define(varMaskName,varMask)
+            #print("+++ defined new var + mask combinations:",varMaskName,varMask)
+            #print(rdf.GetDefinedColumnNames())
+            #print("+++")
+        return rdf,varMaskName
+
+class RDFWrapper:
+    ''' Wrapper of RDataFrame with possibility to define variable+mask combinations
+    '''
+    def __init__(self,rdf):
+        #
+        # variable+mask name indexed by variable+selection string
+        #
+        self.rdf = rdf
+        self.varMaskDefinitions = { }
+
+    def __call__(self):
+        ''' Easy access to RDataFrame
+        '''
+        return self.rdf
+
+    def defineVarMask(self,varName,selection):
+        ''' Return the name of a variable+mask combination. Define if necessary.
+            Arguments:
+              rdf ........ RDataFrame to be used for definition
+              varName .... name of the variable
+              selection .. selection string to be used as mask
+        '''
+        #
+        # use variable name and normalized selection string (remove spaces)
+        #
+        selNorm = selection.replace(" ","")
+        varMask = "(" + varName + ")["+selNorm+"]"
+        if varMask in self.varMaskDefinitions:
+            # combination exists (assume that it's also defined in the RDataFrame)
+            varMaskName = self.varMaskDefinitions[varMask]
+            assert varMaskName in self.rdf.GetDefinedColumnNames()
+        else:
+            # create name for combination and define it in the RDataFrame
+            n = len(self.varMaskDefinitions)
+            varMaskName = "varMask{:03d}".format(n)
+            self.varMaskDefinitions[varMask] = varMaskName
+            self.rdf = self.rdf.Define(varMaskName,varMask)
+            #print("+++ defined new var + mask combinations:",varMaskName,varMask)
+            #print(rdf.GetDefinedColumnNames())
+            #print("+++")
+        return varMaskName
+
+
+def loadConfiguration(configName,selectedNames=[],vetoedNames=[]):
+    ''' Load variable and histogram definitions from a configuration file. The configuration file
+        is expected to contain a "variables" and a "histograms" section. The "histogram" section
         defines dictionaries with the definitions as variable. The name of the variable
         serves as name of the HistogramDefinition.
         Arguments:
            configName ..... name of the module to import from / python file name
            selectedNames .. explicit list of histogram names to be imported (filename wildcard syntax)
            vetoedNames .... explicit list of histogram names to be skipped (filename wildcard syntax)
+        Result:
+           Tuple with dictionary of variable definitions and HistogramDefinitions object
     '''
     # load histogram definitions
     #
-    result = HistogramDefinitions()
+    vDefs = { }
+    hDefs = HistogramDefinitions()
     if configName==None:
-        return result
+        return ( vDefs, hDefs )
 
-    moduleName = configName[:-3] if configName.endswith(".py") else configName
-    module = __import__(moduleName)
-    for n in dir(module):
-        if n.startswith('__'):
-            continue
-        hDict = getattr(module,n)
-        #print(n,type(hDict))
-        #sys.exit()
+    with open(configName,"rt") as yamlFile:
+        allDicts = yaml.load(yamlFile,Loader=yaml.Loader)
+        yamlFile.close()
+    #
+    # store variable definitions as read from configuration
+    #
+    if 'variables' in allDicts:
+        vDefs = allDicts['variables']
+    #
+    # backward compatibility for histograms: treat full input as "histograms" section
+    #   in case neither variables nor histograms sections are defined
+    #
+    if 'histograms' in allDicts:
+        histDicts = allDicts['histograms']
+    else:
+        assert not ( 'variables' in allDicts )
+        histDicts = allDicts
+    #
+    for n,hDict in histDicts.items():
+        #
         assert type(hDict)==dict
         #
         # check if in list of histograms to be displayed
@@ -199,8 +312,8 @@ def loadHistogramDefinitions(configName,selectedNames=[],vetoedNames=[]):
         # add histogram
         #
         hDef = HistogramDefinition(n,hDict)
-        result.add(hDef)
+        hDefs.add(hDef)
         print("Added",hDef.getParameter('canvasName'))
-
-    return result
+        
+    return ( vDefs, hDefs )