Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-10-25 09:37:46

0001 from __future__ import print_function
0002 from __future__ import absolute_import
0003 #
0004 # Misc functions to manipulate Ecal records
0005 # author: Stefano Argiro
0006 # id: $Id: EcalCondTools.py,v 1.14 2013/05/06 08:29:44 argiro Exp $
0007 #
0008 #
0009 # WARNING: we assume that the list of iovs for a given tag
0010 #          contains one element only, the case of several elements
0011 #          will need to be addressed
0012 
0013 #from pluginCondDBPyInterface import *
0014 from CondCore.Utilities import iovInspector as inspect
0015 from ROOT import TCanvas,TH1F, TH2F, gStyle, TChain, TTree, TLegend, TFile
0016 from . import EcalPyUtils
0017 import sys
0018 from math import sqrt
0019 
0020 def listTags(db):
0021     '''List all available tags for a given db '''
0022     try:
0023         db.startReadOnlyTransaction()
0024         tags = db.allTags()
0025         db.commitTransaction()
0026         for tag in tags.split():
0027             print(tag)
0028     except:
0029         return ""
0030 
0031 def listIovs(db,tag):
0032     '''List all available iovs for a given tag'''
0033 
0034     try :
0035        iov = inspect.Iov(db,tag)
0036        iovlist = iov.list()
0037        print("Available iovs for tag: ",tag)
0038        for p in iovlist:
0039            print("  Since " , p[1], " Till " , p[2])
0040      
0041     except Exception as er :
0042         print(" listIovs exception ",er) 
0043 
0044 def dumpXML(db,tag,since,filename='dump.xml'):
0045     '''Dump record in XML format for a given tag '''
0046     try :
0047        iov = inspect.Iov(db,tag)
0048        db.startReadOnlyTransaction()
0049        Plug = __import__(str(db.payloadModules(tag)[0]))
0050        payload = Plug.Object(db)
0051        listOfIovElem= [iovElem for iovElem in db.iov(tag).elements]
0052        inst = 0
0053        for elem in db.iov(tag).elements :
0054            inst = inst + 1
0055            if str(elem.since())==str(since):
0056                found = inst
0057                break
0058        db.commitTransaction()
0059        payload = inspect.PayLoad(db, tag, elem)
0060        out = open(filename,'w')
0061        print(payload, file=out)
0062       
0063     except Exception as er :
0064         print(" dumpXML exception ",er)
0065 
0066 def plot (db, tag,since,filename='plot.root'):
0067     '''Invoke the plot function from the wrapper and save to the specified \
0068        file. The file format will reflect the extension given.'''
0069    
0070     try :
0071        iov = inspect.Iov(db,tag)
0072        db.startReadOnlyTransaction()
0073        Plug = __import__(str(db.payloadModules(tag)[0]))
0074        payload = Plug.Object(db)
0075        listOfIovElem= [iovElem for iovElem in db.iov(tag).elements]
0076        inst = 0
0077        for elem in db.iov(tag).elements :
0078            inst = inst + 1
0079            if str(elem.since())==str(since):
0080                found = inst
0081                break
0082        db.commitTransaction()
0083        payload = inspect.PayLoad(db, tag, elem)
0084        payload.plot(filename,"",[],[])
0085             
0086     except Exception as er :
0087         print(" plot exception ",er)
0088         
0089 
0090 def compare(tag1,db1,since1,
0091             tag2,db2,since2,filename='compare.root'):
0092   '''Produce comparison plots for two records. Save plots to file \
0093      according to format. tag can be an xml file'''
0094   print("EcalCondTools.py compare tag1 ", tag1, "since1 : ", since1, " tag2 ", tag2," s2 ", since2)
0095 
0096   coeff_1_b=[]
0097   coeff_2_b=[]
0098 
0099   coeff_1_e=[]
0100   coeff_2_e=[]
0101   
0102   if  tag1.find(".xml") < 0:
0103       found=0
0104       try:
0105         db1.startReadOnlyTransaction()
0106         Plug = __import__(str(db1.payloadModules(tag1)[0]))
0107         payload = Plug.Object(db1)
0108         listOfIovElem= [iovElem for iovElem in db1.iov(tag1).elements]
0109         what = {'how':'barrel'}
0110         w = inspect.setWhat(Plug.What(),what)
0111         exb = Plug.Extractor(w)
0112         what = {'how':'endcap'}
0113         w = inspect.setWhat(Plug.What(),what)
0114         exe = Plug.Extractor(w)
0115 #        p = getObject(db1,tag1,since1,ex)
0116 #        p.extract(ex)
0117         print(" before loop")
0118         for elem in db1.iov(tag1).elements :       
0119           if str(elem.since())==str(since1):
0120             found=1
0121             payload.load(elem)
0122             payload.extract(exb)
0123             coeff_1_b = [i for i in exb.values()]# first set of coefficients
0124             payload.extract(exe)
0125             coeff_1_e = [i for i in exe.values()]
0126         db1.commitTransaction()
0127 
0128       except Exception as er :
0129         print(" compare first set exception ",er)
0130       if not found :
0131         print("Could not retrieve payload for tag: " , tag1, " since: ", since1)
0132         sys.exit(0)
0133 
0134   else:
0135       coeff_1_b,coeff_1_e = EcalPyUtils.fromXML(tag1)
0136 
0137   if  tag2.find(".xml")<0:
0138       found=0
0139       try:  
0140         db2.startReadOnlyTransaction()
0141         Plug = __import__(str(db2.payloadModules(tag2)[0]))
0142         what = {'how':'barrel'}
0143         w = inspect.setWhat(Plug.What(),what)
0144         exb = Plug.Extractor(w)
0145         what = {'how':'endcap'}
0146         w = inspect.setWhat(Plug.What(),what)
0147         exe = Plug.Extractor(w)
0148 #        p = getObject(db2,tag2,since2)
0149 #        p.extract(ex)
0150         for elem in db2.iov(tag2).elements :       
0151           if str(elem.since())==str(since2):
0152             found=1
0153             payload.load(elem)
0154             payload.extract(exb)
0155             coeff_2_b = [i for i in exb.values()]# second set of coefficients
0156             payload.extract(exe)
0157             coeff_2_e = [i for i in exe.values()]
0158         db2.commitTransaction()
0159      
0160       except Exception as er :
0161           print(" compare second set exception ",er)
0162       if not found :
0163         print("Could not retrieve payload for tag: " , tag2, " since: ", since2)
0164         sys.exit(0)
0165 
0166   else:
0167       coeff_2_b,coeff_2_e = EcalPyUtils.fromXML(tag2)
0168 
0169   gStyle.SetPalette(1)    
0170 
0171   savefile = TFile(filename,"RECREATE")
0172 
0173   ebhisto,ebmap, profx, profy= compareBarrel(coeff_1_b,coeff_2_b)
0174   eephisto,eepmap,eemhisto,eemmap=compareEndcap(coeff_1_e,coeff_2_e)
0175 
0176 #make more canvas (suppressed : use a root file)
0177 
0178 #  cEBdiff = TCanvas("EBdiff","EBdiff")
0179 #  cEBdiff.Divide(2,2)
0180 
0181 #  cEBdiff.cd(1)
0182 #  ebhisto.Draw()
0183 #  cEBdiff.cd(2)
0184 #  ebmap.Draw("colz")
0185 #  cEBdiff.cd(3)
0186 #  profx.Draw()
0187 #  cEBdiff.cd(4)
0188 #  profy.Draw()
0189 
0190 #  EBfilename = "EB_"+filename
0191 #  cEBdiff.SaveAs(EBfilename)
0192 
0193 #  cEEdiff = TCanvas("EEdiff","EEdiff")
0194 #  cEEdiff.Divide(2,2)
0195   
0196 #  cEEdiff.cd(1)
0197 #  eephisto.Draw()
0198 #  cEEdiff.cd(2)
0199 #  eepmap.Draw("colz")
0200 #  cEEdiff.cd(3)
0201 #  eemhisto.Draw()
0202 #  cEEdiff.cd(4)
0203 #  eemmap.Draw("colz")
0204 
0205 #  EEfilename = "EE_"+filename
0206 #  cEEdiff.SaveAs(EEfilename)
0207 
0208 
0209   eeborderphisto,eeborderpmap,eebordermhisto,eebordermmap=compareEndcapBorder(coeff_1_e,coeff_2_e)
0210   ebborderhisto,ebbordermap = compareBarrelBorder(coeff_1_b,coeff_2_b)
0211 
0212   border_diff_distro_can = TCanvas("border_difference","borders difference")
0213   border_diff_distro_can.Divide(2,3)
0214 
0215   border_diff_distro_can.cd(1)
0216   ebborderhisto.Draw()
0217   border_diff_distro_can.cd(2)
0218   ebbordermap.Draw("colz")
0219   border_diff_distro_can.cd(3)
0220   eeborderphisto.Draw()
0221   border_diff_distro_can.cd(4)
0222   eeborderpmap.Draw("colz")
0223   border_diff_distro_can.cd(5)
0224   eebordermhisto.Draw()
0225   border_diff_distro_can.cd(6)
0226   eebordermmap.Draw("colz")
0227 
0228   bordersfilename = "borders_"+filename
0229   prof_filename = "profiles_"+filename
0230   
0231 #  border_diff_distro_can.SaveAs(bordersfilename)
0232 
0233   savefile.Write()
0234 
0235 
0236 
0237 def histo (db, tag,since,filename='histo.root'):
0238     '''Make histograms and save to file. tag can be an xml file'''
0239     
0240     coeff_barl=[]
0241     coeff_endc=[]
0242     
0243     if  tag.find(".xml")< 0:
0244       found=0
0245       try:  
0246 #          exec('import '+db.moduleName(tag)+' as Plug')
0247         db.startReadOnlyTransaction()
0248         Plug = __import__(str(db.payloadModules(tag)[0]))
0249         payload = Plug.Object(db)
0250         listOfIovElem= [iovElem for iovElem in db.iov(tag).elements]
0251         what = {'how':'barrel'}
0252         w = inspect.setWhat(Plug.What(),what)
0253         exb = Plug.Extractor(w)
0254         what = {'how':'endcap'}
0255         w = inspect.setWhat(Plug.What(),what)
0256         exe = Plug.Extractor(w)
0257         for elem in db.iov(tag).elements :       
0258           if str(elem.since())==str(since):
0259             found=1
0260             payload.load(elem)
0261             payload.extract(exb)
0262             coeff_barl = [i for i in exb.values()]
0263             payload.extract(exe)
0264             coeff_endc = [i for i in exe.values()]
0265         db.commitTransaction()
0266 
0267       except Exception as er :
0268           print(" histo exception ",er)
0269       if not found :
0270         print("Could not retrieve payload for tag: " , tag, " since: ", since)
0271         sys.exit(0)
0272 
0273     else :
0274       coeff_barl,coeff_endc=EcalPyUtils.fromXML(tag)
0275 
0276     savefile = TFile(filename,"RECREATE")
0277 
0278     ebmap, ebeta, ebphi, eePmap, ebdist, eeMmap, prof_eePL, prof_eePR, prof_eeML, prof_eeMR, ebBorderdist = makedist(coeff_barl, coeff_endc)
0279 
0280     gStyle.SetPalette(1)
0281 
0282     c =  TCanvas("CCdist")
0283     c.Divide(2,2)
0284 
0285     c.cd(1)  
0286     ebmap.Draw("colz")
0287     c.cd(2)
0288     eePmap.Draw("colz")
0289     c.cd(3)
0290     ebdist.Draw()
0291     c.cd(4)
0292     eeMmap.Draw("colz")
0293 
0294 #    c.SaveAs(filename)
0295 
0296     prof_eb_eta = ebeta.ProfileX()
0297     prof_eb_phi = ebphi.ProfileX()
0298 
0299     c2 = TCanvas("CCprofiles")
0300     c2.Divide(2,2)
0301 
0302     leg = TLegend(0.1,0.7,0.48,0.9)
0303 
0304     c2.cd(1)
0305     prof_eb_eta.Draw()
0306     c2.cd(2)
0307     prof_eb_phi.Draw()
0308     c2.cd(3)
0309     prof_eePL.SetMarkerStyle(8)
0310     prof_eePL.Draw()
0311     prof_eePR.SetMarkerStyle(8)
0312     prof_eePR.SetMarkerColor(2)
0313     prof_eePR.Draw("same")
0314     prof_eeML.SetMarkerStyle(8)
0315     prof_eeML.SetMarkerColor(3)
0316     prof_eeML.Draw("same")
0317     prof_eeMR.SetMarkerStyle(8)
0318     prof_eeMR.SetMarkerColor(5)
0319     prof_eeMR.Draw("same")
0320     leg.AddEntry(prof_eePL,"Dee1","lp")
0321     leg.AddEntry(prof_eePR,"Dee2","lp")
0322     leg.AddEntry(prof_eeMR,"Dee3","lp")
0323     leg.AddEntry(prof_eeML,"Dee4","lp")
0324     leg.Draw()
0325     c2.cd(4)
0326     ebBorderdist.Draw()
0327 
0328     extrafilename = "profiles_"+filename
0329  #   c2.SaveAs(extrafilename)
0330 
0331     savefile.Write()      
0332 
0333 
0334 def getToken(db,tag,since):
0335     ''' Return payload token for a given iov, tag, db'''
0336     try :
0337        iov = inspect.Iov(db,tag)
0338        iovlist = iov.list()
0339        for p in iovlist:
0340            tmpsince=p[1]
0341            if str(tmpsince)==str(since) :
0342                return p[0]
0343        print("Could not retrieve token for tag: " , tag, " since: ", since)
0344        sys.exit(0)
0345        
0346     except Exception as er :
0347        print(er)
0348 
0349 
0350 def getObject(db,tag,since):
0351     ''' Return payload object for a given iov, tag, db'''
0352     found=0
0353     try:
0354 #       exec('import '+db.moduleName(tag)+' as Plug')  
0355        db.startReadOnlyTransaction()
0356        Plug = __import__(str(db.payloadModules(tag)[0]))
0357        print(" getObject Plug")
0358        payload = Plug.Object(db)
0359        db.commitTransaction()
0360        listOfIovElem= [iovElem for iovElem in db.iov(tag).elements]
0361        print(" getObject before loop")
0362        for elem in db.iov(tag).elements :       
0363            if str(elem.since())==str(since):
0364                found=1
0365                print(" getObject found ", elem.since())
0366 #               return Plug.Object(elem)
0367                return elem
0368            
0369     except Exception as er :
0370         print(" getObject exception ",er)
0371 
0372     if not found :
0373         print("Could not retrieve payload for tag: " , tag, " since: ", since)
0374         sys.exit(0)
0375 
0376 
0377 def makedist(coeff_barl, coeff_endc) :
0378 
0379     ebmap = TH2F("EB","EB",360,1,261,171, -85,86)
0380     eePmap = TH2F("EE","EE",100, 1,101,100,1,101)
0381     eeMmap = TH2F("EEminus","EEminus",100,1,101,100,1,101)
0382     ebdist = TH1F("EBdist","EBdist",100,-2,2)
0383     ebBorderdist = TH1F("EBBorderdist","EBBorderdist",100,-2,2)
0384 
0385     ebeta = TH2F("ebeta","ebeta",171,-85,86,100,-2,2)
0386     ebphi = TH2F("ebphi","ebphi",360,1,361,100,-2,2)
0387 
0388     eePL = TH2F("EEPL","EEPlus Left",50,10,55,100,-2,2)
0389     eePR = TH2F("EEPR","EEPlus Right",50,10,55,100,-2,2)
0390     eeML = TH2F("EEML","EEMinus Left",50,10,55,100,-2,2)
0391     eeMR = TH2F("EEMR","EEMinus Right",50,10,55,100,-2,2)
0392     
0393     for i,c in enumerate(coeff_barl):
0394         ieta,iphi = EcalPyUtils.unhashEBIndex(i)
0395         ebmap.Fill(iphi,ieta,c)
0396         ebdist.Fill(c)
0397         ebeta.Fill(ieta,c)
0398         ebphi.Fill(iphi,c)
0399 
0400         if (abs(ieta)==85 or abs(ieta)==65 or abs(ieta)==64 or abs(ieta)==45 or abs(ieta)==44 or abs(ieta)==25 or abs(ieta)==24 or abs(ieta)==1 or iphi%20==1 or iphi%20==0):
0401             ebBorderdist.Fill(c)
0402 
0403 
0404     for i,c in enumerate(coeff_endc):
0405         ix,iy,iz = EcalPyUtils.unhashEEIndex(i)
0406         R = sqrt((ix-50)*(ix-50)+(iy-50)*(iy-50))
0407 
0408         if  iz>0:
0409             eePmap.Fill(ix,iy,c)
0410             if ix<50:
0411                 eePL.Fill(R,c,1)
0412             if ix>50:
0413                 eePR.Fill(R,c,1)
0414 
0415         if iz<0:
0416             eeMmap.Fill(ix,iy,c)
0417             if ix<50:
0418                 eeML.Fill(R,c,1)
0419             if ix>50:
0420                 eeMR.Fill(R,c,1)
0421 
0422     prof_eePL = eePL.ProfileX()
0423     prof_eePR = eePR.ProfileX()
0424     prof_eeML = eeML.ProfileX()
0425     prof_eeMR = eeMR.ProfileX()
0426     
0427     return ebmap, ebeta, ebphi, eePmap, ebdist, eeMmap, prof_eePL, prof_eePR, prof_eeML, prof_eeMR, ebBorderdist
0428 
0429 def compareBarrel(coeff_barl_1,coeff_barl_2) :
0430   '''Return an histogram and a map of the differences '''
0431 
0432   diff_distro_h   = TH1F("diffh","diffh",100,-2,2)
0433   diff_distro_m   = TH2F("diffm","diffm",360,1,361,171,-85,86)
0434   diff_distro_m.SetStats(0)
0435   ebeta = TH2F("ebeta","ebeta",171,-85,86,100,-2,2)
0436   ebphi = TH2F("ebphi","ebphi",360,1,361,100,-2,2)
0437 
0438   
0439   for i,c in enumerate(coeff_barl_1):  
0440       diff = c - coeff_barl_2[i]      
0441       ieta,iphi= EcalPyUtils.unhashEBIndex(i)
0442       diff_distro_h.Fill(diff) 
0443       diff_distro_m.Fill(iphi,ieta,diff)
0444       ebeta.Fill(ieta,diff)
0445       ebphi.Fill(iphi,diff)
0446 
0447   prof_x_h = ebeta.ProfileX()
0448   prof_y_h = ebphi.ProfileX()
0449           
0450   return diff_distro_h, diff_distro_m, prof_x_h, prof_y_h
0451 
0452 
0453 
0454 def compareBarrelBorder(coeff_barl_1,coeff_barl_2) :
0455   '''Return an histogram and a map of the differences '''
0456 
0457   diff_distro_border_h   = TH1F("diffborderh","diffh",100,-2,2)
0458   diff_distro_border_m   = TH2F("diffborderm","diffm",360,1,361,171,-85,86)
0459   diff_distro_border_m.SetStats(0)
0460   
0461   for i,c in enumerate(coeff_barl_1):  
0462       diff = c - coeff_barl_2[i]      
0463       ieta,iphi= EcalPyUtils.unhashEBIndex(i)
0464       if (abs(ieta)==85 or abs(ieta)==65 or abs(ieta)==64 or abs(ieta)==45 or abs(ieta)==44 or abs(ieta)==25 or abs(ieta)==24 or abs(ieta)==1 or iphi%20==1 or iphi%20==0):
0465           diff_distro_border_h.Fill(diff) 
0466       if (abs(ieta)==85 or abs(ieta)==65 or abs(ieta)==64 or abs(ieta)==45 or abs(ieta)==44 or abs(ieta)==25 or abs(ieta)==24 or abs(ieta)==1 or iphi%20==0 or iphi%20==1): 
0467           diff_distro_border_m.Fill(iphi,ieta,diff)
0468           
0469   return diff_distro_border_h, diff_distro_border_m 
0470 
0471 
0472 
0473 
0474     
0475 def compareEndcap(coeff_endc_1, coeff_endc_2) :
0476     ''' Return an histogram and a map of the differences for each endcap'''
0477 
0478     diff_distro_h_eep   = TH1F("diff EE+","diff EE+",100,-2,2)
0479     diff_distro_h_eem   = TH1F("diff EE-","diff EE-",100,-2,2)
0480 
0481     
0482     diff_distro_m_eep   = TH2F("map EE+","map EE+",100,1,101,100,1,101)
0483     diff_distro_m_eem   = TH2F("map EE-","map EE-",100,1,101,100,1,101)
0484 
0485     temp_h = TH1F("tempR","tempR",50,0,50)
0486     
0487     diff_distro_m_eep.SetStats(0)
0488     diff_distro_m_eem.SetStats(0)
0489 
0490 
0491     for i,c in enumerate(coeff_endc_1):  
0492       diff = c - coeff_endc_2[i]
0493       ix,iy,iz = EcalPyUtils.unhashEEIndex(i)
0494       R = sqrt((ix-50)*(ix-50)+(iy-50)*(iy-50))
0495       
0496       if iz >0:
0497           diff_distro_h_eep.Fill(diff)
0498           diff_distro_m_eep.Fill(ix,iy,diff)
0499           
0500       else:
0501           diff_distro_h_eem.Fill(diff)
0502           diff_distro_m_eem.Fill(ix,iy,diff)
0503 
0504     return diff_distro_h_eep, \
0505            diff_distro_m_eep, \
0506            diff_distro_h_eem, \
0507            diff_distro_m_eem
0508 
0509 
0510 
0511 def compareEndcapBorder(coeff_endc_1, coeff_endc_2) :
0512     ''' Return an histogram and a map of the differences for each endcap'''
0513 
0514     border_diff_distro_h_eep   = TH1F("borderdiff EE+","diff EE+",100,-2,2)
0515     border_diff_distro_h_eem   = TH1F("borderdiff EE-","diff EE-",100,-2,2)
0516 
0517     
0518     border_diff_distro_m_eep   = TH2F("bordermap EE+","map EE+",100,1,101,100,1,101)
0519     border_diff_distro_m_eem   = TH2F("bordermap EE-","map EE-",100,1,101,100,1,101)
0520     
0521     border_diff_distro_m_eep.SetStats(0)
0522     border_diff_distro_m_eem.SetStats(0)
0523 
0524 
0525     for i,c in enumerate(coeff_endc_1):  
0526       diff = c - coeff_endc_2[i]
0527       ix,iy,iz = EcalPyUtils.unhashEEIndex(i)
0528       Rsq = ((ix-50.0)**2+(iy-50.0)**2)
0529       
0530       if (iz >0 and (Rsq<144.0 or Rsq>2500.0)):
0531           border_diff_distro_h_eep.Fill(diff)
0532           border_diff_distro_m_eep.Fill(ix,iy,diff)
0533       elif (iz<0 and (Rsq<144.0 or Rsq>2500.0)):
0534           border_diff_distro_h_eem.Fill(diff)
0535           border_diff_distro_m_eem.Fill(ix,iy,diff)
0536       
0537 
0538     return border_diff_distro_h_eep, \
0539            border_diff_distro_m_eep, \
0540            border_diff_distro_h_eem, \
0541            border_diff_distro_m_eem