| # Gaudi Framework, version v21r8 |

Home |  | Generated: 17 Mar 2010

# compareRootHistos.py

Go to the documentation of this file.

    00001 from ROOT import TFile
    00002 import sets
    00003 import sys
    00004
    00005 histos = ['TH1D', 'TH2D', 'TProfile']
    00006 ser = 'SERIAL'
    00007 par = 'PARALL'
    00008
    00009 # =================================================================================================
    00010 # Method   : rec( o, path=None, lst=None )
    00011 #
    00012 # @param o    : a ROOT object
    00013 # @param path : a string like a transient store path; ie '/stat/CaloPIDs/ECALPIDE'
    00014 # @param lst  : a list to hold (path, object) tuples
    00015 #
    00016 # function    : recursively pull apart a ROOT file, making a list of (path, TObject) tuples
    00017 #               This is done by GetListOfKeys method, which lets one work down through directories
    00018 #               until you hit the Histo at the end of the path.  The list of tuples is returned
    00019 #
    00020 def rec( o, path=None, lst=None ) :
    00021   if not path : path = '/stat' ; lst = []
    00022   else        : path = path + '/' + o.GetName()
    00023   lst.append( (path,o) )
    00024   if 'GetListOfKeys' in dir(o) :
    00025     keys = o.GetListOfKeys()
    00026     for k in keys :
    00027       name =  k.GetName()
    00028       rec( o.Get(name), path, lst )
    00029   else :
    00030     pass
    00031   return lst
    00032 # =================================================================================================
    00033
    00034 # =================================================================================================
    00035 # Method   : composition( t )
    00036 #
    00037 # @param t : a tuple of ( type, d ) where type is either 'SERIAL' or 'PARALL'
    00038 #            and d is a dictionary of ROOT objects, with each key = ROOT path
    00039 #
    00040 # function : deduce the composition, (objects/histos) counts
    00041 #
    00042 def composition( t ) :
    00043   typ, d = t
    00044   hists = 0 ; objs = 0
    00045   for k in d.keys() :
    00046     if d[k].__class__.__name__ in histos : hists += 1
    00047     else                                 : objs  += 1
    00048   return objs, hists
    00049 # =================================================================================================
    00050
    00051 # =================================================================================================
    00052 # Method        : comparePaths( t1, t2 )
    00053 #
    00054 # @param t1, t2 : a tuple of ( type, d ) where type is either 'SERIAL' or 'PARALL'
    00055 #                 and d is a dictionary of ROOT objects, with each key = ROOT path
    00056 #
    00057 # function      : compare the paths between the two histo files.  If the files are identical, they
    00058 #                 should have the same set of paths.  The Parallel file should definitely have the
    00059 #                 same paths as the Serial.  Perhaps the Serial file will have some more paths due
    00060 #                 to extra histos added as part of Application Sequencer finalisation
    00061 #                 Arguments t1 and t2 are checked and the parallel/serial auto-detected
    00062 #                 Uses sets module for intersections/unions, etc.
    00063 #
    00064 def comparePaths( t1, t2 ) :
    00065   if   t1[0]  == ser : ds = t1[1] ; dp = t2[1]
    00066   elif t2[0]  == ser : ds = t2[1] ; dp = t1[1]
    00067   else : print 'Neither tuple is Serial Root file reference?' ; return
    00068
    00069   dsks = ds.keys() ; dpks = dp.keys()
    00070   dsks.sort()      ; dpks.sort()
    00071
    00072   sset = sets.Set( dsks )
    00073   pset = sets.Set( dpks )
    00074   os, hs = composition( (ser, ds) )
    00075   op, hp = composition( (par, dp) )
    00076   print '\n' + '='*80
    00077   print 'Comparison of Paths : Serial vs Parallel ROOT files'
    00078   print '-'*80
    00079   print 'Number of paths in Serial file : %i (objects, histos) = ( %i, %i )'%( len(dsks), os, hs )
    00080   print 'Number of paths in Parall file : %i (objects, histos) = ( %i, %i )'%( len(dpks), op, hp )
    00081   matching = sset.intersection(pset)
    00082   matchingHistos = 0
    00083   for n in matching :
    00084     if ds[n].__class__.__name__ in histos : matchingHistos += 1
    00085   print '\nMatching paths                 : %i'%( len(matching) )
    00086   uSer = sset - pset
    00087   # work out histos unique to parallel file
    00088   uniqueSerialHistos = 0
    00089   for n in uSer :
    00090     if ds[n].__class__.__name__ in histos : uniqueSerialHistos += 1
    00091   print 'Paths unique to Serial file : %i ( %i Histos )'%( len(uSer), uniqueSerialHistos )
    00092   if uSer :
    00093     for n in uSer : print '\t%s : \t%s'%( ds[n], n )
    00094   uPar = pset - sset
    00095   uniqueParallHistos = 0
    00096   for n in uPar :
    00097     if dp[n].__class__.__name__ in histos : uniqueParallHistos += 1
    00098   print 'Paths unique to Parall file : %i ( %i Histos )'%( len(uPar), uniqueParallHistos )
    00099   if uPar :
    00100     for n in uPar : print '\t%s : \t%s'%( dp[n], n )
    00101   print 'Matching Histos to test : %i'%( matchingHistos )
    00102   print '='*80 + '\n'
    00103   return ( ((os,hs),(op,hp)), (uSer, uniqueSerialHistos), (uPar, uniqueParallHistos), matchingHistos )
    00104 # =================================================================================================
    00105
    00106 # =================================================================================================
    00107 # Method        : compareHistos( t1, t2 )
    00108 #
    00109 # @param t1, t2 : a tuple of ( type, d ) where type is either 'SERIAL' or 'PARALL'
    00110 #                 and d is a dictionary of ROOT objects, with each key = ROOT path
    00111 #
    00112 # function      : compare the histograms in Serial/Parallel ROOT files.  First, go through each
    00113 #                 dict to collect the histos (ignore TDirectory objects, etc).  Then the histos
    00114 #                 in the parallel file (experimental) are compared to their equivalents in the
    00115 #                 serial file (definitely correct) using 3 methods.
    00116 #                 1) The entries are checked, they should be equal
    00117 #                 2) If entries are equal, check the Integral(); should be equal
    00118 #                 3) If integrals are equal, check the KolmogorovTest() ; should be 1
    00119 #                 Arguments t1 and t2 are checked and the parallel/serial auto-detected
    00120 #
    00121 def compareHistos(t1, t2, state) :
    00122
    00123   ( ((serialObjects,serialHistos),(parallObjects, parallHistos)), (uniqueSerPaths,uniqueSerHistos), (uniqueParPaths,uniqueParHistos), mh ) = state
    00124
    00125   # deduce which one is parallel, which serial
    00126   if   t1[0]  == ser : ds = t1[1] ; dp = t2[1]
    00127   elif t2[0]  == ser : ds = t2[1] ; dp = t1[1]
    00128   else : print 'Neither tuple is Serial Root file reference?' ; return
    00129
    00130   # histocount, objectcount for parallel/serial
    00131   hcp = 0 ; pHistos = []
    00132   hcs = 0 ; sHistos = []
    00133
    00134   omit = ['/stat/Brunel/MemoryTool/Virtual mem, all entries',
    00135           '/stat/Brunel/MemoryTool/Virtual mem, downscaled']
    00136   omit = []
    00137
    00138   # find the histos in the serial file
    00139   for k in ds.keys() :
    00140     if k not in omit :
    00141       if ds[k].__class__.__name__ in histos : hcs += 1 ; sHistos.append( k )
    00142   # same for parallel
    00143   for k in dp.keys() :
    00144     if k not in omit :
    00145       if dp[k].__class__.__name__ in histos : hcp += 1 ; pHistos.append( k )
    00146
    00147
    00148   cEntries   = 0 ; xEntries   = 0 ; diffEntries   = []
    00149   cIntegrals = 0 ; xIntegrals = 0 ; diffIntegrals = []
    00150   passedKol  = 0 ; failedKol  = 0 ; diffKols      = [] ; zeroIntegrals = 0
    00151   kTested    = 0
    00152   notfound   = 0 ; integralMatch = 0 ; otherTest = 0 ; zeroIntegralMatch = 0
    00153   for h in sHistos :
    00154     if h in pHistos :
    00155       # matching histos to check
    00156       cEntries += 1
    00157       sh = ds[h] ; ph = dp[h]
    00158       # first check entries
    00159       if sh.GetEntries() != ph.GetEntries() : diffEntries.append(h) ; xEntries += 1 ; continue
    00160       # check for (non-zero sum of bin error) && (non-zero integrals) for K-Test
    00161       sBinError = 0.0 ; pBinError = 0.0
    00162       for i in xrange(sh.GetNbinsX()) : sBinError += sh.GetBinError(i)
    00163       for i in xrange(ph.GetNbinsX()) : pBinError += ph.GetBinError(i)
    00164       sint = sh.Integral() ; pint = ph.Integral()
    00165       if (bool(sint) and bool(pint)) and ( sBinError>0 and pBinError>0 ) :
    00166         kTested += 1
    00167         kTest = sh.KolmogorovTest(ph)
    00168         if int(kTest) : passedKol += 1
    00169         else          : failedKol += 1 ; diffKols.append(h) # ; print 'KTest result : ', kTest
    00170       else :
    00171         # try the integral test?
    00172         otherTest += 1
    00173         if all((sint, pint)) and (sint==pint) :
    00174           integralMatch += 1
    00175         elif (sint==pint) :
    00176           zeroIntegralMatch += 1
    00177         else :
    00178             diffIntegrals.append( h )
    00179             xIntegrals += 1
    00180     else :
    00181       notfound += 1 ; print 'not found? ', h
    00182
    00183   # report on Failed Entry-Checks
    00184   print '\n\n'+'-'*80
    00185   print 'Summary of histos with different Entries'
    00186   print '-'*80
    00187   if diffEntries :
    00188     diffEntries.sort()
    00189     for e in diffEntries : print '\t\t\t%s:\t%i != %i'%( e, int(ds[e].GetEntries()), int(dp[e].GetEntries()) )
    00190   print '-'*80
    00191
    00192   # report on Failed Kolmogorov Tests
    00193   print '\n\n'+'-'*60
    00194   print 'Summary of histos which failed Kolmogorov Test'
    00195   print '-'*60
    00196   if diffKols :
    00197     diffKols.sort()
    00198     for e in diffKols :
    00199       result = ds[e].KolmogorovTest(dp[e])
    00200       print '%s\t\t%s :\tK-Test Result :\t %5.16f'%( type(ds[e]), e, result )
    00201   print '-'*60
    00202
    00203   # report on Failed Integral Checks
    00204   print '\n\n'+'-'*60
    00205   print 'Summary of histos which failed Integral Check'
    00206   print '-'*60
    00207   if diffIntegrals :
    00208     diffIntegrals.sort()
    00209     for e in diffIntegrals :
    00210       diff = dp[e].Integral()-ds[e].Integral()
    00211       pc   = (diff*100)/ds[e].Integral()
    00212       print '%s\t\t%s:\t Diff = %5.6f\tPercent Diff to Serial : %5.6f '%( type(ds[e]), e, diff, pc  )
    00213   print '-'*60 + '\n'
    00214   print '='*80 + '\n'
    00215
    00216   print '\n' + '='*80
    00217   print 'Comparison : Serial/Parallel ROOT Histo files'
    00218   print '\n\t\tSerial\tParall'
    00219   print '\tObjects : %i\t%i\t\t( p-s = %i )'%( serialObjects, parallObjects, parallObjects-serialObjects )
    00220   print '\tHistos  : %i\t%i\t\t( p-s = %i )'%( serialHistos,  parallHistos,  parallHistos-serialHistos  )
    00221   print '\t          __________'
    00222   print '\tTotal   : %i\t%i\n'%( serialHistos+serialObjects,  parallHistos+parallObjects  )
    00223   print 'Objects/Histos unique to Serial File : %i / %i'%( len(uniqueSerPaths)-uniqueSerHistos, uniqueSerHistos )
    00224   print 'Objects/Histos unique to Parall File : %i / %i'%( len(uniqueParPaths)-uniqueParHistos, uniqueParHistos )
    00225   print '\nMatching Histograms valid for Comparison : %i'%( mh )
    00226   print '\nOmissions : '
    00227   for entry in omit : print '\t%s'%( entry )
    00228   print '\nHistograms for Comparison (after Omissions) : %i'%( mh-len(omit) )
    00229   print '\n\tHISTOGRAM TESTS : '
    00230   print '\t\tKOLMOGOROV TEST      : %i'%( kTested )
    00231   print '\t\tINTEGRAL TEST        : %i'%( otherTest )
    00232   print '\t\tENTRIES TEST         : %i'%( xEntries )
    00233   print '\t\t                       ____'
    00234   print '\t\tTested               : %i'%( cEntries )
    00235
    00236   print '\n\tDISCREPANCIES : '
    00237   print '\t\tK-Test      : %i'%( failedKol  )
    00238   print '\t\tIntegrals   : %i'%( xIntegrals )
    00239   print '\t\tEntries     : %i'%( xEntries   )
    00240   print '\n'+'='*80
    00241
    00242 # =================================================================================================
    00243
    00244 if __name__ == '__main__' :
    00245   sys.argv.pop(0)   # get rid of script name
    00246   if len(sys.argv) == 2 :
    00247     pFile = sys.argv[0]
    00248     sFile = sys.argv[1]
    00249   else :
    00250     print '*'*80
    00251     print 'Wrong count of arguments? > python compareRootHistos.py someParallelFile.root someSerialFile.root'
    00252     print '*'*80
    00253     sys.exit(0)
    00254   tfs = TFile( sFile, 'REC' ) ; print 'opening Serial File : %s'%( sFile )
    00255   tfp = TFile( pFile, 'REC' ) ; print 'opening Parall File : %s'%( pFile )
    00256
    00257   # get structure of TFiles in a list of (path, object) tuples
    00258   lser = rec(tfs) ; lpar = rec(tfp)
    00259   # make a dictionary of lser and lpar.  keys=paths
    00260   dserial = dict( [(n, o) for n, o in lser] )
    00261   dparall = dict( [(n, o) for n, o in lpar] )
    00262   # make a tuple of (type, dict) where type is either 'serial' or 'parallel'
    00263   ts = ( ser, dserial )  ; tp = ( par, dparall )
    00264
    00265   # check objs/histos in each file
    00266 #   composition( ts ) ; composition( tp )
    00267
    00268   # compare paths from each file
    00269   state = comparePaths( ts, tp )
    00270
    00271   # compare histos from each file
    00272   compareHistos( ts, tp, state )
    00273
    00274 #   # finished with TFiles
    00275 #   tfs.Close()     ; tfp.Close()
    00276 #
    00277