Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-11-28 23:11:15

0001 #! /usr/bin/env python3
0002 
0003 # Pure trick to start ROOT in batch mode, pass this only option to it
0004 # and the rest of the command line options to this code.
0005 import sys
0006 import numpy as np
0007 import pandas as pd
0008 import matplotlib.pyplot as plt
0009 
0010 from array import array
0011 oldargv = sys.argv[:]
0012 sys.argv = [ '-b-' ]
0013 from ROOT import TCanvas, TLegend, TPaveText, THStack, TFile, TLatex, TPaveLabel
0014 from ROOT import TProfile, TProfile2D, TH1D, TH2D, TH2F, TPaletteAxis, TH1, TH1F, TColor, TExec, TLine
0015 from ROOT import kBlack, kWhite, kOrange, kAzure, kBlue
0016 from ROOT import gROOT, gStyle
0017 gROOT.SetBatch(True)
0018 sys.argv = oldargv
0019 
0020 from Validation.Geometry.plot_utils import setTDRStyle, Plot_params, drawEtaValues
0021  
0022 from Validation.Geometry.plot_hgcal_utils import plots, COMPOUNDS, DETECTORS, sDETS, hist_label_to_num, TwikiPrintout, acustompalette, dEdx, MatXo, drawHalfEtaValues
0023 from collections import namedtuple, OrderedDict
0024 import sys, os
0025 import argparse
0026 
0027 def paramsGood_(detector, plot):
0028     """Check the validity of the arguments.
0029 
0030        Common function to check the validity of the parameters passed
0031        in. It returns a tuple composed by a bool and a string. The
0032        bool indicates if all checks are ok, the string the name of the
0033        appropriate ROOT file to open (empty string in case the any
0034        check failed)
0035 
0036     """
0037 
0038     if plot not in plots.keys():
0039         print("Error, unknown plot %s" % plot)
0040         return (False, '')
0041 
0042     if detector not in DETECTORS and detector not in COMPOUNDS.keys():
0043         print('Error, unknown detector: %s' % detector)
0044         return (False, '')
0045 
0046     theDetectorFilename = ''
0047     if detector in DETECTORS:
0048         theDetectorFilename = 'matbdg_%s.root' % detector
0049     else:
0050         theDetectorFilename = 'matbdg_%s.root' % COMPOUNDS[detector][0]
0051 
0052     if not checkFile_(theDetectorFilename):
0053         print("Error, missing file %s" % theDetectorFilename)
0054         raise RuntimeError
0055     return (True, theDetectorFilename)
0056 
0057 def checkFile_(filename):
0058     return os.path.exists(filename)
0059 
0060 def setColorIfExists_(histos, h, color):
0061     if h in histos.keys():
0062         histos[h].SetFillColor(color)
0063 
0064 def assignOrAddIfExists_(h, p):
0065     """Assign the projection of p to h.
0066 
0067        Function to assign the projection of p to h, in the case in
0068        which h is None, otherwise add the projection to the already
0069        valid h object
0070 
0071     """
0072 
0073     if not h:
0074         h = p.ProjectionX()
0075     else:
0076         h.Add(p.ProjectionX("B_%s" % h.GetName()), +1.000)
0077     return h
0078 
0079 def createPlots_(plot, compounddetectorname):
0080     """Cumulative material budget from simulation.
0081     
0082        Internal function that will produce a cumulative profile of the
0083        material budget inferred from the simulation starting from the
0084        single detectors that compose the tracker. It will iterate over
0085        all existing detectors contained in the DETECTORS
0086        dictionary. The function will automatically skip non-existent
0087        detectors.
0088 
0089     """
0090 
0091     theDirname = "Figures"
0092 
0093     hist_X0_detectors = OrderedDict()
0094     if plot not in plots.keys():
0095         print("Error: chosen plot name not known %s" % plot)
0096         return
0097 
0098     # We need to keep the file content alive for the lifetime of the
0099     # full function....
0100     subDetectorFiles = []
0101 
0102     hist_X0_elements = OrderedDict()
0103     prof_X0_elements = OrderedDict()
0104     for subDetector,color in DETECTORS.items():
0105         subDetectorFilename = "matbdg_%s.root" % subDetector
0106         if not checkFile_(subDetectorFilename):
0107             print("Error opening file: %s" % subDetectorFilename)
0108             continue
0109 
0110         subDetectorFiles.append(TFile(subDetectorFilename))
0111         subDetectorFile = subDetectorFiles[-1]
0112         print ("Opening file: %s" % subDetectorFilename)
0113         prof_X0_XXX = subDetectorFile.Get("%d" % plots[plot].plotNumber)
0114 
0115         hist_X0_detectors[subDetector] = prof_X0_XXX.ProjectionX()
0116 
0117         # category profiles
0118         for label, [num, color, leg] in hist_label_to_num.items():
0119             prof_X0_elements[label] = subDetectorFile.Get("%d" % (num + plots[plot].plotNumber))
0120             hist_X0_elements[label] = assignOrAddIfExists_(hist_X0_elements.setdefault(label, None),
0121                                                           prof_X0_elements[label])
0122 
0123     cumulative_matbdg = TH1D("CumulativeSimulMatBdg",
0124                              "CumulativeSimulMatBdg",
0125                              hist_X0_detectors["BeamPipe"].GetNbinsX(),
0126                              hist_X0_detectors["BeamPipe"].GetXaxis().GetXmin(),
0127                              hist_X0_detectors["BeamPipe"].GetXaxis().GetXmax())
0128     cumulative_matbdg.SetDirectory(0)
0129 
0130     # colors
0131     for det, color in DETECTORS.items():
0132         setColorIfExists_(hist_X0_detectors, det, color)
0133 
0134     for label, [num, color, leg] in hist_label_to_num.items():
0135         hist_X0_elements[label].SetFillColor(color)
0136 
0137     # First Plot: BeamPipe + Tracker + ECAL + HCal + HGCal + MB + MGNT
0138     # stack
0139     stackTitle_SubDetectors = "Material Budget;%s;%s" % (
0140         plots[plot].abscissa,plots[plot].ordinate)
0141     stack_X0_SubDetectors = THStack("stack_X0",stackTitle_SubDetectors)
0142     for det, histo in hist_X0_detectors.items():
0143         stack_X0_SubDetectors.Add(histo)
0144         cumulative_matbdg.Add(histo, 1)
0145 
0146     # canvas
0147     can_SubDetectors = TCanvas("can_SubDetectors","can_SubDetectors",800,800)
0148     #can_SubDetectors.Range(0,0,25,25)
0149     can_SubDetectors.SetFillColor(kWhite)
0150 
0151     # Draw
0152     stack_X0_SubDetectors.SetMinimum(plots[plot].ymin)
0153     stack_X0_SubDetectors.SetMaximum(plots[plot].ymax)
0154     stack_X0_SubDetectors.Draw("HIST")
0155     #stack_X0_SubDetectors.GetXaxis().SetLimits(plots[plot].xmin, plots[plot].xmax)
0156 
0157 
0158     # Legenda
0159     theLegend_SubDetectors = TLegend(0.130,0.7,0.93,0.90) #(0.180,0.8,0.98,0.90)
0160     theLegend_SubDetectors.SetNColumns(2)
0161     theLegend_SubDetectors.SetFillColor(0)
0162     theLegend_SubDetectors.SetFillStyle(0)
0163     theLegend_SubDetectors.SetBorderSize(0)
0164 
0165     for det, histo in hist_X0_detectors.items():
0166         theLegend_SubDetectors.AddEntry(histo, det,  "f")
0167 
0168     theLegend_SubDetectors.Draw()
0169 
0170     # text
0171     text_SubDetectors = TPaveText(0.130,0.627,0.352,0.687,"NDC") #(0.180,0.727,0.402,0.787,"NDC")
0172     text_SubDetectors.SetFillColor(0)
0173     text_SubDetectors.SetBorderSize(0)
0174     text_SubDetectors.AddText("CMS Simulation")
0175     text_SubDetectors.SetTextAlign(11)
0176     text_SubDetectors.Draw()
0177 
0178     # Store
0179     can_SubDetectors.Update()
0180     if not checkFile_(theDirname):
0181         os.mkdir(theDirname)
0182     can_SubDetectors.SaveAs("%s/MaterialBdg_%s_%s.pdf" % (theDirname, compounddetectorname,plot))
0183     can_SubDetectors.SaveAs("%s/MaterialBdg_%s_%s.png" % (theDirname, compounddetectorname,plot))
0184     can_SubDetectors.SaveAs("%s/MaterialBdg_%s_%s.root" % (theDirname, compounddetectorname,plot))
0185 
0186     if plot == "x_vs_eta" or plot == "l_vs_eta":
0187         canname = "MBCan_1D_%s_%s_total"  % (compounddetectorname, plot)
0188         can2 = TCanvas(canname, canname, 800, 800)
0189         can2.Range(0,0,25,25)
0190         can2.SetFillColor(kWhite)
0191         gStyle.SetOptStat(0)
0192         gStyle.SetOptTitle(0);
0193         #title = TPaveLabel(.11,.95,.35,.99,"Total accumulated material budget","brndc")
0194         stack_X0_SubDetectors.GetStack().Last().SetMarkerStyle(34)
0195         stack_X0_SubDetectors.GetStack().Last().GetXaxis().SetRangeUser( 1.0, 3.5)
0196         stack_X0_SubDetectors.GetStack().Last().Draw();
0197         stack_X0_SubDetectors.GetYaxis().SetTitleOffset(1.15);
0198         can2.Update()
0199         can2.Modified()
0200         can2.SaveAs("%s/%s_%s_total_Zplus.pdf" % (theDirname, compounddetectorname, plot))
0201         can2.SaveAs("%s/%s_%s_total_Zplus.png" % (theDirname, compounddetectorname, plot))
0202         stack_X0_SubDetectors.GetStack().Last().GetXaxis().SetRangeUser( -3.5, -1.0)
0203         stack_X0_SubDetectors.GetStack().Last().Draw();
0204         stack_X0_SubDetectors.GetYaxis().SetTitleOffset(1.15);
0205         can2.Update()
0206         can2.Modified()
0207         can2.SaveAs("%s/%s_%s_total_Zminus.pdf" % (theDirname, compounddetectorname, plot))
0208         can2.SaveAs("%s/%s_%s_total_Zminus.png" % (theDirname, compounddetectorname, plot))
0209 
0210         #Also print them to give them exact numbers
0211         etavalues = []
0212         matbudginX0 = []
0213         matbudginIntLen = []
0214         for binx in range(0,stack_X0_SubDetectors.GetStack().Last().GetXaxis().GetNbins()):
0215             bincontent = stack_X0_SubDetectors.GetStack().Last().GetBinContent(binx) 
0216             if bincontent == 0: continue
0217             etavalues.append( stack_X0_SubDetectors.GetStack().Last().GetBinCenter(binx)  )
0218             if plot == "x_vs_eta": 
0219                 matbudginX0.append( bincontent  )
0220                 d1 = {'Eta': etavalues, 'MatBudInX0': matbudginX0}
0221                 df1 = pd.DataFrame(data=d1).round(2)
0222                 df1.to_csv(r'/afs/cern.ch/work/a/apsallid/CMS/PFCalStudies/CMS-HGCAL/matbudV10fromVertexToBackofHGCal/CMSSW_11_0_X_2019-06-04-2300/src/Validation/Geometry/test/EtavsMatBudinXo.txt',sep=' ', index=False, header=False)
0223                 #print df1
0224             if plot == "l_vs_eta": 
0225                 matbudginIntLen.append( bincontent )
0226                 d2 = {'Eta': etavalues, 'MatBudInIntLen': matbudginIntLen}
0227                 df2 = pd.DataFrame(data=d2).round(2)
0228                 df2.to_csv(r'/afs/cern.ch/work/a/apsallid/CMS/PFCalStudies/CMS-HGCAL/matbudV10fromVertexToBackofHGCal/CMSSW_11_0_X_2019-06-04-2300/src/Validation/Geometry/test/EtavsMatBudInIntLen.txt',sep=' ', index=False, header=False)
0229                 #print df2
0230 
0231     return cumulative_matbdg
0232 
0233 def createPlots2D_(plot, compounddetectorname):
0234     """2D material budget map to know exactly what we are adding.
0235     """
0236 
0237     #IBs = ["InnerServices", "Phase2PixelBarrel", "TIB", "TIDF", "TIDB"]
0238     theDirname = "Figures"
0239 
0240     hist_X0_detectors = OrderedDict()
0241     if plot not in plots.keys():
0242         print("Error: chosen plot name not known %s" % plot)
0243         return
0244 
0245     # We need to keep the file content alive for the lifetime of the
0246     # full function....
0247     subDetectorFiles = []
0248 
0249     hist_X0_elements = OrderedDict()
0250     prof_X0_elements = OrderedDict()
0251 
0252     for subDetector,color in DETECTORS.items():
0253         subDetectorFilename = "matbdg_%s.root" % subDetector
0254         if not checkFile_(subDetectorFilename):
0255             print("Error opening file: %s" % subDetectorFilename)
0256             continue
0257 
0258         subDetectorFiles.append(TFile(subDetectorFilename))
0259         subDetectorFile = subDetectorFiles[-1]
0260         print ("Opening file: %s" % subDetectorFilename)
0261         prof_X0_XXX = subDetectorFile.Get("%d" % plots[plot].plotNumber)
0262 
0263         #hist_X0_detectors[subDetector] = prof_X0_XXX
0264         hist_X0_detectors[subDetector] = prof_X0_XXX.ProjectionXY("_pxy","B")
0265         print(subDetector)
0266 
0267     # First Plot: BeamPipe + Tracker + ECAL + HCal + HGCal + MB + MGNT
0268  
0269     # Create "null" histo
0270     minX = 1.03*hist_X0_detectors["BeamPipe"].GetXaxis().GetXmin()
0271     maxX = 1.03*hist_X0_detectors["BeamPipe"].GetXaxis().GetXmax()
0272     minY = 1.03*hist_X0_detectors["BeamPipe"].GetYaxis().GetXmin()
0273     maxY = 1.03*hist_X0_detectors["BeamPipe"].GetYaxis().GetXmax()
0274 
0275     frame = TH2F("frame", "", 10, minX, maxX, 10, minY, maxY);
0276     frame.SetMinimum(0.1)
0277     frame.SetMaximum(10.)
0278     frame.GetXaxis().SetTickLength(frame.GetXaxis().GetTickLength()*0.50)
0279     frame.GetYaxis().SetTickLength(frame.GetXaxis().GetTickLength()/4.)
0280 
0281     hist2d_X0_total = hist_X0_detectors["BeamPipe"]
0282 
0283     # stack
0284     hist2dTitle = ('%s %s;%s;%s;%s' % (plots[plot].quotaName,
0285                                        "All detectors",
0286                                        plots[plot].abscissa,
0287                                        plots[plot].ordinate,
0288                                        plots[plot].quotaName))
0289 
0290     hist2d_X0_total.SetTitle(hist2dTitle)
0291     frame.SetTitle(hist2dTitle)
0292     frame.SetTitleOffset(0.5,"Y")
0293 
0294     #If here you put different histomin,histomaxin plot_utils you won't see anything 
0295     #for the material plots. 
0296     if plots[plot].histoMin != -1.:
0297         hist2d_X0_total.SetMinimum(plots[plot].histoMin)
0298     if plots[plot].histoMax != -1.:
0299         hist2d_X0_total.SetMaximum(plots[plot].histoMax)
0300 
0301     #
0302     # canvas
0303     can_SubDetectors = TCanvas("can_SubDetectors","can_SubDetectors",2480+248, 580+58+58)
0304     can_SubDetectors.SetTopMargin(0.1)
0305     can_SubDetectors.SetBottomMargin(0.1)
0306     can_SubDetectors.SetLeftMargin(0.04)
0307     can_SubDetectors.SetRightMargin(0.06)
0308     can_SubDetectors.SetFillColor(kWhite)
0309     gStyle.SetOptStat(0)
0310     gStyle.SetTitleFillColor(0)
0311     gStyle.SetTitleBorderSize(0)
0312     gStyle.SetOptTitle(0)
0313 
0314     hist2d_X0_total.GetYaxis().SetTickLength(hist2d_X0_total.GetXaxis().GetTickLength()/4.)
0315     hist2d_X0_total.GetYaxis().SetTickLength(hist2d_X0_total.GetXaxis().GetTickLength()/4.)
0316     hist2d_X0_total.SetTitleOffset(0.5,"Y")
0317     hist2d_X0_total.GetYaxis().SetTitleOffset(0.50);
0318     #hist2d_X0_total.GetXaxis().SetTitleOffset(1.15);
0319     #hist2d_X0_total.GetXaxis().SetNoExponent(True)
0320     #hist2d_X0_total.GetYaxis().SetNoExponent(True)
0321 
0322 
0323     # colors
0324     for det, color in DETECTORS.items():
0325         hist_X0_detectors[det].SetMarkerColor(color)
0326         hist_X0_detectors[det].SetFillColor(color)
0327 
0328     for det, histo in hist_X0_detectors.items():
0329         print(det)
0330         histo.Draw("same")
0331 
0332     # Legenda
0333     theLegend_SubDetectors = TLegend(0.100,0.7,0.90,0.90)#(0.180,0.8,0.98,0.90)
0334     theLegend_SubDetectors.SetNColumns(3)
0335     theLegend_SubDetectors.SetFillColor(0)
0336     theLegend_SubDetectors.SetFillStyle(0)
0337     theLegend_SubDetectors.SetBorderSize(0)
0338 
0339     for det, histo in hist_X0_detectors.items():
0340         theLegend_SubDetectors.AddEntry(histo, det,  "f")
0341     #theLegend_SubDetectors.AddEntry(hgbound1, "HGCal Eta Boundaries [1.3, 3.0]",  "l") 
0342 
0343     theLegend_SubDetectors.Draw()
0344 
0345     
0346     # text
0347     text_SubDetectors = TPaveText(0.100,0.627,0.322,0.687,"NDC")#(0.180,0.727,0.402,0.787,"NDC")
0348     text_SubDetectors.SetFillColor(0)
0349     text_SubDetectors.SetBorderSize(0)
0350     text_SubDetectors.AddText("CMS Simulation")
0351     text_SubDetectors.SetTextAlign(11)
0352     text_SubDetectors.Draw()
0353     
0354 
0355     #Add eta labels
0356     keep_alive = []
0357     if plots[plot].iDrawEta:
0358         keep_alive.extend(drawEtaValues())
0359     
0360     # Store
0361     can_SubDetectors.Update()
0362     if not checkFile_(theDirname):
0363         os.mkdir(theDirname)
0364     can_SubDetectors.SaveAs("%s/MaterialBdg_%s_%s.png" % (theDirname, compounddetectorname, plot))
0365     #It seems that it is too heavy to create .pdf and .root
0366     #can_SubDetectors.SaveAs("%s/MaterialBdg_FromVertexToEndofHGCal_%s.pdf" % (theDirname, plot))
0367     #can_SubDetectors.SaveAs("%s/MaterialBdg_FromVertexToEndofHGCal_%s.root" % (theDirname, plot))
0368 
0369     hist2d_X0_total.GetXaxis().SetRangeUser( 0., 7000.)
0370     #Draw eta values in the zoom case
0371     keep_alive = []
0372     keep_alive.extend(drawHalfEtaValues())
0373     #hist2d_X0_total.Draw("COLZ") 
0374     can_SubDetectors.Update()
0375     can_SubDetectors.Modified()
0376     can_SubDetectors.SaveAs("%s/MaterialBdg_%s_%s_Zpluszoom.png" % (theDirname, compounddetectorname, plot))
0377 
0378 def createPlotsReco_(reco_file, label, debug=False):
0379     """Cumulative material budget from reconstruction.
0380 
0381 
0382        Internal function that will produce a cumulative profile of the
0383        material budget in the reconstruction starting from the single
0384        detectors that compose the tracker. It will iterate over all
0385        existing detectors contained in the sDETS dictionary. The
0386        function will automatically stop everytime it encounters a
0387        non-existent detector, until no more detectors are left to
0388        try. For this reason the keys in the sDETS dictionary can be as
0389        inclusive as possible.
0390 
0391     """
0392 
0393     cumulative_matbdg = None
0394     sPREF = ["Original_RadLen_vs_Eta_", "RadLen_vs_Eta_"]
0395 
0396     c = TCanvas("c", "c", 1024, 1024);
0397     diffs = []
0398     if not checkFile_(reco_file):
0399         print("Error: missing file %s" % reco_file)
0400         raise RuntimeError
0401     file = TFile(reco_file)
0402     prefix = "/DQMData/Run 1/Tracking/Run summary/RecoMaterial/"
0403     for s in sPREF:
0404         hs = THStack("hs","");
0405         histos = []
0406         for det, color in sDETS.items():
0407             layer_number = 0
0408             while True:
0409                 layer_number += 1
0410                 name = "%s%s%s%d" % (prefix, s, det, layer_number)
0411                 prof = file.Get(name)
0412                 # If we miss an object, since we are incrementally
0413                 # searching for consecutive layers, we may safely
0414                 # assume that there are no additional layers and skip
0415                 # to the next detector.
0416                 if not prof:
0417                     if debug:
0418                         print("Missing profile %s" % name)
0419                     break
0420                 else:
0421                     histos.append(prof.ProjectionX("_px", "hist"));
0422                     diffs.append(histos[-1]);
0423                     histos[-1].SetFillColor(color + layer_number);
0424                     histos[-1].SetLineColor(color + layer_number + 1);
0425 
0426         name = "CumulativeRecoMatBdg_%s" % s
0427         if s == "RadLen_vs_Eta_":
0428             cumulative_matbdg = TH1D(name, name,
0429                                      histos[0].GetNbinsX(),
0430                                      histos[0].GetXaxis().GetXmin(),
0431                                      histos[0].GetXaxis().GetXmax())
0432             cumulative_matbdg.SetDirectory(0)
0433         for h in histos:
0434             hs.Add(h)
0435             if cumulative_matbdg:
0436                 cumulative_matbdg.Add(h, 1.)
0437         hs.Draw()
0438         hs.GetYaxis().SetTitle("RadLen")
0439         c.Update()
0440         c.Modified()
0441         c.SaveAs("%sstacked_%s.png" % (s, label))
0442     hs = THStack("diff","")
0443     for d in range(0,len(diffs)/2):
0444         diffs[d+len(diffs)/2].Add(diffs[d], -1.)
0445         hs.Add(diffs[d+len(diffs)/2]);
0446     hs.Draw()
0447     hs.GetYaxis().SetTitle("RadLen")
0448     c.Update()
0449     c.Modified()
0450     c.SaveAs("RadLen_difference_%s.png" % label)
0451     return cumulative_matbdg
0452 
0453 def materialBudget_Simul_vs_Reco(reco_file, label, debug=False):
0454     """Plot reco vs simulation material budget.
0455     
0456        Function are produces a direct comparison of the material
0457        budget as extracted from the reconstruction geometry and
0458        inferred from the simulation one.
0459 
0460     """
0461 
0462     setTDRStyle()
0463 
0464     # plots
0465     cumulative_matbdg_sim = createPlots_("x_vs_eta")
0466     cumulative_matbdg_rec = createPlotsReco_(reco_file, label, debug=False)
0467 
0468     cc = TCanvas("cc", "cc", 1024, 1024)
0469     cumulative_matbdg_sim.SetMinimum(0.)
0470     cumulative_matbdg_sim.SetMaximum(3.5)
0471     cumulative_matbdg_sim.GetXaxis().SetRangeUser(-3.0, 3.0)
0472     cumulative_matbdg_sim.SetLineColor(kOrange)
0473     cumulative_matbdg_rec.SetMinimum(0.)
0474     cumulative_matbdg_rec.SetMaximum(3.)
0475     cumulative_matbdg_rec.SetLineColor(kAzure+1)
0476     l = TLegend(0.18, 0.8, 0.95, 0.92)
0477     l.AddEntry(cumulative_matbdg_sim, "Sim Material", "f")
0478     l.AddEntry(cumulative_matbdg_rec, "Reco Material", "f")
0479     cumulative_matbdg_sim.Draw("HIST")
0480     cumulative_matbdg_rec.Draw("HIST SAME")
0481     l.Draw()
0482     filename = "MaterialBdg_Reco_vs_Simul_%s.png" % label
0483     cc.SaveAs(filename)
0484 
0485 def createCompoundPlots(detector, plot):
0486     """Produce the requested plot for the specified detector.
0487 
0488        Function that will plot the requested @plot for the specified
0489        @detector. The specified detector could either be a real
0490        detector or a compound one. The list of available plots are the
0491        keys of plots dictionary (imported from plot_utils.
0492 
0493     """
0494 
0495     theDirname = 'Images'
0496     if not checkFile_(theDirname):
0497         os.mkdir(theDirname)
0498 
0499     goodToGo, theDetectorFilename = paramsGood_(detector, plot)
0500     if not goodToGo:
0501         return
0502 
0503     theDetectorFile = TFile(theDetectorFilename)
0504     #
0505 
0506     # get TProfiles
0507     prof_X0_elements = OrderedDict()
0508     hist_X0_elements = OrderedDict()
0509     for label, [num, color, leg] in hist_label_to_num.items():
0510         #print label, num, color, leg
0511         prof_X0_elements[label] = theDetectorFile.Get("%d" % (num + plots[plot].plotNumber))
0512         hist_X0_elements[label] = prof_X0_elements[label].ProjectionX()
0513         hist_X0_elements[label].SetFillColor(color)
0514         hist_X0_elements[label].SetLineColor(kBlack)
0515 
0516     files = []
0517     if detector in COMPOUNDS.keys():
0518         for subDetector in COMPOUNDS[detector][1:]:
0519             subDetectorFilename = "matbdg_%s.root" % subDetector
0520 
0521             # open file
0522             if not checkFile_(subDetectorFilename):
0523                 continue
0524 
0525             subDetectorFile = TFile(subDetectorFilename)
0526             files.append(subDetectorFile)
0527             print("*** Open file... %s" %  subDetectorFilename)
0528 
0529             # subdetector profiles
0530             for label, [num, color, leg] in hist_label_to_num.items():
0531                 prof_X0_elements[label] = subDetectorFile.Get("%d" % (num + plots[plot].plotNumber))
0532                 hist_X0_elements[label].Add(prof_X0_elements[label].ProjectionX("B_%s" % prof_X0_elements[label].GetName())
0533                                             , +1.000)
0534 
0535     # stack
0536     stackTitle = "Material Budget %s;%s;%s" % (detector,
0537                                                plots[plot].abscissa,
0538                                                plots[plot].ordinate)
0539     stack_X0 = THStack("stack_X0", stackTitle);
0540     for label, [num, color, leg] in hist_label_to_num.items():
0541         stack_X0.Add(hist_X0_elements[label])
0542 
0543     # canvas
0544     canname = "MBCan_1D_%s_%s"  % (detector, plot)
0545     can = TCanvas(canname, canname, 800, 800)
0546     can.Range(0,0,25,25)
0547     can.SetFillColor(kWhite)
0548     gStyle.SetOptStat(0)
0549     gStyle.SetOptTitle(1);
0550 
0551     # Draw
0552     stack_X0.Draw("HIST");
0553     stack_X0.GetYaxis().SetTitleOffset(1.15);
0554 
0555     # Legenda
0556     theLegend = TLegend(0.40, 0.65, 0.60, 0.89)
0557     if plot == "x_vs_phi" or plot == "l_vs_phi": theLegend = TLegend(0.65, 0.30, 0.89, 0.70)
0558     if plot == "x_vs_R" or plot == "l_vs_R": theLegend = TLegend(0.75, 0.60, 0.95, 0.90)
0559 
0560     for label, [num, color, leg] in hist_label_to_num.items():
0561         theLegend.AddEntry(hist_X0_elements[label], leg, "f")
0562     theLegend.Draw();
0563 
0564     # Store
0565     can.Update();
0566     can.SaveAs( "%s/%s_%s.pdf" % (theDirname, detector, plot))
0567     can.SaveAs( "%s/%s_%s.png" % (theDirname, detector, plot))
0568 
0569     #Let's also save the total accumulated budget vs eta since muon id relies 
0570     #on adequate calorimeter thickness
0571     if plot == "x_vs_eta" or plot == "l_vs_eta":
0572         canname = "MBCan_1D_%s_%s_total"  % (detector, plot)
0573         can2 = TCanvas(canname, canname, 800, 800)
0574         can2.Range(0,0,25,25)
0575         can2.SetFillColor(kWhite)
0576         gStyle.SetOptStat(0)
0577         gStyle.SetOptTitle(0);
0578         #title = TPaveLabel(.11,.95,.35,.99,"Total accumulated material budget","brndc")
0579         stack_X0.GetStack().Last().GetXaxis().SetRangeUser( 0., 3.)
0580         stack_X0.GetStack().Last().Draw();
0581         stack_X0.GetYaxis().SetTitleOffset(1.15);
0582         can2.Update()
0583         can2.Modified()
0584         can2.SaveAs("%s/%s_%s_total_Zplus.pdf" % (theDirname, detector, plot))
0585         can2.SaveAs("%s/%s_%s_total_Zplus.png" % (theDirname, detector, plot))
0586         stack_X0.GetStack().Last().GetXaxis().SetRangeUser( -3., 0.)
0587         stack_X0.GetStack().Last().Draw();
0588         stack_X0.GetYaxis().SetTitleOffset(1.15);
0589         can2.Update()
0590         can2.Modified()
0591         can2.SaveAs("%s/%s_%s_total_Zminus.pdf" % (theDirname, detector, plot))
0592         can2.SaveAs("%s/%s_%s_total_Zminus.png" % (theDirname, detector, plot))
0593  
0594 
0595 def create2DPlots(detector, plot, plotnum, plotmat, dosingledetector = True):
0596     """Produce the requested plot for the specified detector.
0597 
0598        Function that will plot the requested 2D-@plot for the
0599        specified @detector. The specified detector could either be a
0600        real detector or a compound one. The list of available plots
0601        are the keys of plots dictionary imported from plot_utils.
0602     """
0603 
0604     #gStyle.Reset()
0605     #Better to use an underscore. 
0606     plotmat = plotmat.replace(" ", "_")
0607 
0608 
0609     if plotmat != "": 
0610         theDirname = ('Images/%s' % plotmat).replace(" ", "")
0611     else: 
0612         theDirname = 'Images'
0613 
0614     if not checkFile_(theDirname):
0615         os.mkdir(theDirname)
0616     if not os.path.isdir(('Images/%s/ZPlusZoom' % plotmat).replace(" ", "")):
0617         os.mkdir( ('Images/%s/ZPlusZoom' % plotmat).replace(" ", "") )
0618     if not os.path.isdir(('Images/%s/ZMinusZoom' % plotmat).replace(" ", "")):
0619         os.mkdir( ('Images/%s/ZMinusZoom' % plotmat).replace(" ", "") )
0620 
0621     goodToGo, theDetectorFilename = paramsGood_(detector, plot)
0622     if not goodToGo:
0623         return
0624 
0625     theDetectorFile = TFile(theDetectorFilename)
0626 
0627     prof2d_X0_det_total = TProfile2D()
0628     prof2d_X0_det_total.Reset()
0629 
0630     # get TProfiles
0631     #prof2d_X0_det_total = theDetectorFile.Get('%s' % plots[plot].plotNumber)
0632     prof2d_X0_det_total = theDetectorFile.Get('%s' % plotnum)
0633     print("==================================================================") 
0634     print(plotnum)
0635 
0636     # histos
0637     prof2d_X0_det_total.__class__ = TProfile2D
0638     hist_X0_total = prof2d_X0_det_total.ProjectionXY()
0639 
0640     # keep files live forever
0641     files = []
0642     if detector in COMPOUNDS.keys() and not dosingledetector:
0643         #When the loop was:
0644         #for subDetector in COMPOUNDS[detector][1:]:
0645         #and the detector was single it never went in the loop and read the single file 
0646         #from above. I alter this to COMPOUNDS[detector] to do the multi material budget plot.
0647         #This won't effect the single detector due to the alter in the if above
0648         for subDetector in COMPOUNDS[detector]:
0649             # filenames of single components
0650             subDetectorFilename = "matbdg_%s.root" % subDetector
0651  
0652             # open file
0653             if not checkFile_(subDetectorFilename):
0654                 print("Error, missing file %s" % subDetectorFilename)
0655                 continue
0656     
0657             subDetectorFile = TFile(subDetectorFilename)
0658             files.append(subDetectorFile)
0659             print("*** Open file... %s" %  subDetectorFilename)
0660 
0661             # subdetector profiles
0662             prof2d_X0_det_total = subDetectorFile.Get('%s' % plots[plot].plotNumber)
0663             prof2d_X0_det_total.__class__ = TProfile2D
0664 
0665             # add to summary histogram
0666             hist_X0_total.Add(prof2d_X0_det_total.ProjectionXY("B_%s" % prof2d_X0_det_total.GetName()), +1.000 )
0667 
0668     # # properties
0669     #gStyle.SetPalette(1)
0670     gStyle.SetStripDecimals(False)
0671     # #
0672 
0673     # Create "null" histo
0674     minX = 1.03*prof2d_X0_det_total.GetXaxis().GetXmin()
0675     maxX = 1.03*prof2d_X0_det_total.GetXaxis().GetXmax()
0676     minY = 1.03*prof2d_X0_det_total.GetYaxis().GetXmin()
0677     maxY = 1.03*prof2d_X0_det_total.GetYaxis().GetXmax()
0678 
0679     frame = TH2F("frame", "", 10, minX, maxX, 10, minY, maxY);
0680     frame.SetMinimum(0.1)
0681     frame.SetMaximum(10.)
0682     frame.GetXaxis().SetTickLength(frame.GetXaxis().GetTickLength()*0.50)
0683     frame.GetYaxis().SetTickLength(frame.GetXaxis().GetTickLength()/4.)
0684 
0685     # Ratio
0686     if plots[plot].iRebin:
0687         prof2d_X0_det_total.Rebin2D()
0688 
0689     # stack
0690     hist2dTitle = ('%s %s;%s;%s;%s' % (plots[plot].quotaName,
0691                                        detector,
0692                                        plots[plot].abscissa,
0693                                        plots[plot].ordinate,
0694                                        plots[plot].quotaName))
0695 
0696     if dosingledetector: 
0697         hist2d_X0_total = prof2d_X0_det_total
0698     else : 
0699         hist2d_X0_total = hist_X0_total
0700     hist2d_X0_total.SetTitle(hist2dTitle)
0701     frame.SetTitle(hist2dTitle)
0702     frame.SetTitleOffset(0.5,"Y")
0703 
0704     #If here you put different histomin,histomaxin plot_utils you won't see anything 
0705     #for the material plots. 
0706     if plots[plot].histoMin != -1.:
0707         hist2d_X0_total.SetMinimum(plots[plot].histoMin)
0708     if plots[plot].histoMax != -1.:
0709         hist2d_X0_total.SetMaximum(plots[plot].histoMax)
0710 
0711     #
0712     can2name = "MBCan_2D_%s_%s_%s" % (detector, plot, plotmat)
0713     can2 = TCanvas(can2name, can2name, 2480+248, 580+58+58)
0714     can2.SetTopMargin(0.1)
0715     can2.SetBottomMargin(0.1)
0716     can2.SetLeftMargin(0.04)
0717     can2.SetRightMargin(0.06)
0718     can2.SetFillColor(kWhite)
0719     gStyle.SetOptStat(0)
0720     gStyle.SetTitleFillColor(0)
0721     gStyle.SetTitleBorderSize(0)
0722 
0723     #hist2d_X0_total.SetMaximum(hist2d_X0_total.GetMaximum())
0724     # Color palette
0725     # gStyle.SetPalette()#1
0726     acustompalette()
0727     ex1 = TExec("ex1","acustompalette();");
0728     ex1.Draw();
0729 
0730     #for i in range(100): MyPaletteArray.append(i+1)
0731 
0732     #gStyle.SetPalette(first_color_number);
0733 
0734     # Log?
0735     can2.SetLogz(plots[plot].zLog)
0736 
0737     # Draw in colors
0738     #frame.Draw()
0739     #hist2d_X0_total.Draw("COLZsame") #Dummy draw to create the palette object
0740     hist2d_X0_total.Draw("COLZ") #Dummy draw to create the palette object
0741 
0742     # Store
0743     can2.Update()
0744 
0745     #Aesthetic
0746     palette = hist2d_X0_total.GetListOfFunctions().FindObject("palette")
0747     if palette:
0748         palette.__class__ = TPaletteAxis
0749         palette.SetX1NDC(0.945)
0750         palette.SetX2NDC(0.96)
0751         palette.SetY1NDC(0.1)
0752         palette.SetY2NDC(0.9)
0753         palette.GetAxis().SetTickSize(.01)
0754         palette.GetAxis().SetTitle("")
0755         if plots[plot].zLog:
0756             palette.GetAxis().SetLabelOffset(-0.01)
0757     paletteTitle = TLatex(1.12*maxX, maxY, plots[plot].quotaName)
0758     paletteTitle.SetTextAngle(90.)
0759     paletteTitle.SetTextSize(0.05)
0760     paletteTitle.SetTextAlign(31)
0761     paletteTitle.Draw()
0762     hist2d_X0_total.GetYaxis().SetTickLength(hist2d_X0_total.GetXaxis().GetTickLength()/4.)
0763     hist2d_X0_total.GetYaxis().SetTickLength(hist2d_X0_total.GetXaxis().GetTickLength()/4.)
0764     hist2d_X0_total.SetTitleOffset(0.5,"Y")
0765     hist2d_X0_total.GetYaxis().SetTitleOffset(0.45);
0766     #hist2d_X0_total.GetXaxis().SetTitleOffset(1.15);
0767     #hist2d_X0_total.GetXaxis().SetNoExponent(True)
0768     #hist2d_X0_total.GetYaxis().SetNoExponent(True)
0769 
0770     #Add eta labels
0771     keep_alive = []
0772     if plots[plot].iDrawEta:
0773         keep_alive.extend(drawEtaValues())
0774 
0775     can2.Modified()
0776     hist2d_X0_total.SetContour(255)
0777 
0778     # Store
0779     can2.Update()
0780     can2.Modified()
0781 
0782     can2.SaveAs( "%s/%s_%s%s.pdf" % (theDirname, detector, plot, plotmat))
0783     can2.SaveAs( "%s/%s_%s%s.png" % (theDirname, detector, plot, plotmat))
0784     #can2.SaveAs( "%s/%s_%s%s.root" % (theDirname, detector, plot, plotmat))
0785 
0786     #Zoom in a little bit
0787     if plot == "x_vs_z_vs_Rsum" or plot == "l_vs_z_vs_Rsum" or plot == "x_vs_z_vs_Rsumcos" or plot == "l_vs_z_vs_Rsumcos" or plot == "x_vs_z_vs_Rloc" or plot == "l_vs_z_vs_Rloc" or  plot == "x_vs_z_vs_Rloccos" or plot == "l_vs_z_vs_Rloccos":
0788         #Z+
0789         #hist2d_X0_total.GetXaxis().SetLimits( 3100., 5200.)
0790         if dosingledetector: hist2d_X0_total.GetXaxis().SetRangeUser( 3100., 5400.)
0791         else : hist2d_X0_total.GetXaxis().SetRangeUser( 0., 7000.)
0792         #Do not draw eta values in the zoom case
0793         keep_alive = []
0794         #hist2d_X0_total.Draw("COLZ") 
0795         can2.Update()
0796         can2.Modified()
0797         can2.SaveAs( "%s/%s/%s_%s%s_ZplusZoom.pdf" % (theDirname, "ZPlusZoom", detector, plot, plotmat))
0798         can2.SaveAs( "%s/%s/%s_%s%s_ZplusZoom.png" % (theDirname, "ZPlusZoom", detector, plot, plotmat))
0799         #Z-
0800         #hist2d_X0_total.GetXaxis().SetLimits( 3100., 5200.)
0801         if dosingledetector: hist2d_X0_total.GetXaxis().SetRangeUser( -5400., -3100.)
0802         else : hist2d_X0_total.GetXaxis().SetRangeUser( 0., -7000.)
0803         #Do not draw eta values in the zoom case
0804         keep_alive = []
0805         #hist2d_X0_total.Draw("COLZ") 
0806         can2.Update()
0807         can2.Modified()
0808         can2.SaveAs( "%s/%s/%s_%s%s_ZminusZoom.pdf" % (theDirname, "ZMinusZoom", detector, plot, plotmat))
0809         can2.SaveAs( "%s/%s/%s_%s%s_ZminusZoom.png" % (theDirname, "ZMinusZoom", detector, plot, plotmat))
0810 
0811 
0812     gStyle.SetStripDecimals(True)
0813 
0814 def createRatioPlots(detector, plot):
0815     """Create ratio plots.
0816 
0817        Function that will make the ratio between the radiation length
0818        and interaction length, for the specified detector. The
0819        specified detector could either be a real detector or a
0820        compound one.
0821 
0822     """
0823 
0824     goodToGo, theDetectorFilename = paramsGood_(detector, plot)
0825     if not goodToGo:
0826         return
0827 
0828     theDirname = 'Images'
0829     if not os.path.exists(theDirname):
0830         os.mkdir(theDirname)
0831 
0832     theDetectorFile = TFile(theDetectorFilename)
0833     # get TProfiles
0834     prof_x0_det_total = theDetectorFile.Get('%d' % plots[plot].plotNumber)
0835     prof_l0_det_total = theDetectorFile.Get('%d' % (10000+plots[plot].plotNumber))
0836 
0837     # histos
0838     hist_x0_total = prof_x0_det_total.ProjectionX()
0839     hist_l0_total = prof_l0_det_total.ProjectionX()
0840 
0841     if detector in COMPOUNDS.keys():
0842         for subDetector in COMPOUNDS[detector][1:]:
0843 
0844             # file name
0845             subDetectorFilename = "matbdg_%s.root" % subDetector
0846 
0847             # open file
0848             if not checkFile_(subDetectorFilename):
0849                 print("Error, missing file %s" % subDetectorFilename)
0850                 continue
0851 
0852             subDetectorFile = TFile(subDetectorFilename)
0853 
0854             # subdetector profiles
0855             prof_x0_det_total = subDetectorFile.Get('%d' % plots[plot].plotNumber)
0856             prof_l0_det_total = subDetectorFile.Get('%d' % (10000+plots[plot].plotNumber))
0857             # add to summary histogram
0858             hist_x0_total.Add(prof_x0_det_total.ProjectionX("B_%s" % prof_x0_det_total.GetName()), +1.000 )
0859             hist_l0_total.Add(prof_l0_det_total.ProjectionX("B_%s" % prof_l0_det_total.GetName()), +1.000 )
0860     #
0861     hist_x0_over_l0_total = hist_x0_total
0862     hist_x0_over_l0_total.Divide(hist_l0_total)
0863     histTitle = "Material Budget %s;%s;%s" % (detector,
0864                                               plots[plot].abscissa,
0865                                               plots[plot].ordinate)
0866     hist_x0_over_l0_total.SetTitle(histTitle)
0867     # properties
0868     hist_x0_over_l0_total.SetMarkerStyle(1)
0869     hist_x0_over_l0_total.SetMarkerSize(3)
0870     hist_x0_over_l0_total.SetMarkerColor(kBlue)
0871 
0872     # canvas
0873     canRname = "MBRatio_%s_%s" % (detector, plot)
0874     canR = TCanvas(canRname,canRname,800,800)
0875     canR.Range(0,0,25,25)
0876     canR.SetFillColor(kWhite)
0877     gStyle.SetOptStat(0)
0878 
0879     # Draw
0880     hist_x0_over_l0_total.Draw("E1")
0881 
0882     # Store
0883     canR.Update()
0884     canR.SaveAs("%s/%s_%s.pdf" % (theDirname, detector, plot))
0885     canR.SaveAs("%s/%s_%s.png" % (theDirname, detector, plot))
0886 
0887     
0888 def GetSiliconZValuesFromXML(): 
0889     """Someone can use the xml files of the geometry. Just do 
0890     cd $CMSSW_BASE/src/Geometry/HGCalCommonData/test
0891     cmsRun testHGCalParameters_cfg.py
0892     The scintillator position in the HEB is very close to the silicon but they have a 
0893     1.95 mm difference. So, we take the silicon position. 
0894 
0895     This returns a np.array of the z position of the silicon values in mm. 
0896     """
0897     
0898     #First as they come out of the testHGCalParameters_cfg.py as cm 
0899     #plus either the thermal screen front (297.0 cm)
0900     #or the EE front (319.05 cm) 
0901     layersmaxZ = ( 319.05, 319.815, 320.725, 322.255, 323.165, 324.695, 325.605, 327.135, 328.045, 329.575, 330.485, 332.015, 332.925, 334.455, 335.365, 336.895, 337.805, 339.335, 340.245, 341.775, 342.685, 344.215, 345.125, 346.655, 347.565, 349.095, 350.005, 351.535, 352.445, 357.715, 362.615, 367.515, 372.415, 377.315, 382.215, 387.115, 392.015, 397.105, 402.195, 407.285, 412.375, 420.765, 429.155, 437.545, 445.935, 454.325, 462.715, 471.105, 479.495, 487.885, 496.275, 504.665, 513.055)
0902 
0903     layersmaxZ = np.asarray(layersmaxZ)
0904     layersmaxZ = 10. * layersmaxZ # in mm
0905 
0906     print("Total number of layers from XML " , layersmaxZ.size - 1) # Minus 1 for the EE front input by hand. 
0907 
0908     return layersmaxZ
0909 
0910 
0911 if __name__ == '__main__':
0912     parser = argparse.ArgumentParser(description='Generic Material Plotter',
0913                                      formatter_class=argparse.ArgumentDefaultsHelpFormatter)
0914     parser.add_argument('-r', '--reco',
0915                         help='Input Reconstruction Material file, DQM format')
0916     parser.add_argument('-l', '--label',
0917                         help='Label to use in naming the plots')
0918     parser.add_argument('-c', '--compare',
0919                         help='Compare simulation and reco materials',
0920                         action='store_true',
0921                         default=False)
0922     parser.add_argument('-m', '--multi',
0923                         help='Combining multiple detectors',
0924                         action='store_true',
0925                         default=False)
0926     parser.add_argument('-s', '--single',
0927                         help='Material budget for single detector from simulation',
0928                         action='store_true',
0929                         default=False)
0930     parser.add_argument('-d', '--detector',
0931                         help='Detector for which you want to compute the material budget',
0932                         type=str,)
0933     args = parser.parse_args()
0934 
0935     if args.compare and args.single:
0936         print("Error, too many actions required")
0937         raise RuntimeError
0938 
0939     if args.compare:
0940         if args.reco is None:
0941             print("Error, missing inpur reco file")
0942             raise RuntimeError
0943         if args.label is None:
0944             print("Error, missing label")
0945             raise RuntimeError
0946         materialBudget_Simul_vs_Reco(args.reco, args.label, debug=False)
0947 
0948     if args.multi: 
0949         #Material Budget plot from Vertex to end of HGCal
0950         cumulative_matbdg_sim = createPlots_("x_vs_eta",args.detector)
0951         cumulative_matbdg_sim = createPlots_("l_vs_eta",args.detector)
0952         cumulative_matbdg_sim = createPlots2D_("x_vs_z_vs_Rsum",args.detector)
0953         cumulative_matbdg_sim = createPlots2D_("l_vs_z_vs_Rsum",args.detector)
0954 
0955         required_2DplotsForalldetectors = ["x_vs_z_vs_Rsum", "l_vs_z_vs_Rsum"]
0956 
0957         for p in required_2DplotsForalldetectors:
0958             #First the total
0959             create2DPlots(args.detector, p, plots[p].plotNumber, "", not args.multi)
0960             #Then, the rest
0961             #for label, [num, color, leg] in hist_label_to_num.iteritems():
0962                 #print label, num, color, leg
0963                 #create2DPlots(args.detector, p, num + plots[p].plotNumber, leg)
0964 
0965     if args.single:
0966         if args.detector is None:
0967             print("Error, missing detector")
0968             raise RuntimeError
0969         required_2Dplots = ["x_vs_eta_vs_phi", "l_vs_eta_vs_phi", "x_vs_z_vs_Rsum", "l_vs_z_vs_Rsum", "x_vs_z_vs_Rsumcos", "l_vs_z_vs_Rsumcos", "x_vs_z_vs_Rloc", "l_vs_z_vs_Rloc", "x_vs_z_vs_Rloccos", "l_vs_z_vs_Rloccos"]
0970         # The plots below require z bin 1 cm to have more statistics for eta, phi 
0971         # while the plots above are with 200000 events for z bin 1 mm. 
0972         #required_2Dplots = ["x_vs_eta_vs_phi", "l_vs_eta_vs_phi", "x_vs_z_vs_Rsum", "l_vs_z_vs_Rsum", "x_vs_z_vs_Rsumcos"]
0973 
0974         required_plots = ["x_vs_eta", "x_vs_phi", "x_vs_R", "l_vs_eta", "l_vs_phi", "l_vs_R"]
0975 
0976         required_ratio_plots = ["x_over_l_vs_eta", "x_over_l_vs_phi"]
0977 
0978         #Function to help filling the twiki with all these plots
0979         #First I loop through labels to put the hide button in twiki
0980         #All HGCal
0981         print("---+++ Results: Plots for individual material in all HGCal")
0982         for label, [num, color, leg] in hist_label_to_num.items():
0983             for p in ["x_vs_z_vs_Rsum", "l_vs_z_vs_Rsum", "x_vs_z_vs_Rsumcos", "l_vs_z_vs_Rsumcos", "x_vs_z_vs_Rloc", "l_vs_z_vs_Rloc"]:
0984                 TwikiPrintout(p, leg, "all")
0985         #Z+
0986         print("---+++ Results: Plots for individual material in Z+ Endcap")
0987         for label, [num, color, leg] in hist_label_to_num.items():
0988             for p in ["x_vs_z_vs_Rsum", "l_vs_z_vs_Rsum", "x_vs_z_vs_Rsumcos", "l_vs_z_vs_Rsumcos", "x_vs_z_vs_Rloc", "l_vs_z_vs_Rloc"]:
0989                 TwikiPrintout(p, leg, "zplus")
0990         #Z-
0991         print("---+++ Results: Plots for individual material in Z- Endcap")
0992         for label, [num, color, leg] in hist_label_to_num.items():
0993             for p in ["x_vs_z_vs_Rsum", "l_vs_z_vs_Rsum", "x_vs_z_vs_Rsumcos", "l_vs_z_vs_Rsumcos", "x_vs_z_vs_Rloc", "l_vs_z_vs_Rloc"]:
0994                 TwikiPrintout(p, leg, "zminus")
0995 
0996         #Below is the silicon position from the xml geometry file
0997         #Should we put them on top of plots like the eta values?
0998         print(GetSiliconZValuesFromXML())
0999 
1000         for p in required_2Dplots:
1001             #First the total
1002             create2DPlots(args.detector, p, plots[p].plotNumber, "")
1003             #Then, the rest
1004             for label, [num, color, leg] in hist_label_to_num.items():
1005                 #print label, num, color, leg
1006                 create2DPlots(args.detector, p, num + plots[p].plotNumber, leg)
1007 
1008 
1009         for p in required_plots:
1010             createCompoundPlots(args.detector, p)
1011         for p in required_ratio_plots:
1012             createRatioPlots(args.detector, p)