Ecal/python/EcalCondTools.py

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