Page MenuHomeHEPForge
Files F8308654

gap_removal
No OneTemporary
Actions

gap_removal
View Options

#! /usr/bin/env python
import os, sys, tempfile
class Inputdata:
def __init__(self, filename):
self.histos = {}
self.description = {}
self.description['DrawOnly'] = []
f = open(filename+'.dat', 'r')
for line in f:
if (line.count('#',0,1)):
if (line.count('BEGIN HISTOGRAM')):
title = line.split('BEGIN HISTOGRAM', 1)[1].strip()
self.description['DrawOnly'].append(title)
self.histos[title] = Histogram(f)
f.close()
class Histogram:
def __init__(self, f):
self.read_input(f)
def read_input(self, f):
self.description = {}
self.data = []
for line in f:
if (line.count('#',0,1)):
if (line.count('END HISTOGRAM')):
break
else:
line = line.rstrip()
if (line.count('=')):
linearray = line.split('=', 1)
self.description[linearray[0]] = linearray[1]
else:
linearray = line.split('\t')
if len(linearray)==4:
self.data.append({'LowEdge': float(linearray[0]),
'UpEdge': float(linearray[1]),
'Content': float(linearray[2]),
'Error': [float(linearray[3]),float(linearray[3])]})
elif len(linearray)==5:
self.data.append({'LowEdge': float(linearray[0]),
'UpEdge': float(linearray[1]),
'Content': float(linearray[2]),
'Error': [float(linearray[3]),float(linearray[4])]})
def write_datapoint(self, f, xval, xerr, yval, yerr):
f.write(' <dataPoint>\n')
f.write(' <measurement errorPlus="%e" value="%e" errorMinus="%e"/>\n' %(xerr, xval, xerr))
f.write(' <measurement errorPlus="%e" value="%e" errorMinus="%e"/>\n' %(yerr[1], yval, yerr[0]))
f.write(' </dataPoint>\n')
def write_datapointset_header(self, f):
path = self.description['AidaPath']
title = self.description['Title'].replace('>', '&gt;').replace('<', '&lt;').replace('"', '&quot;')
f.write(' <dataPointSet name="%s" dimension="2"\n' % path.split('/')[-1])
f.write(' path="%s" title="%s">\n' %(os.path.abspath(path.replace(path.split('/')[-1], '')), title))
f.write(' <annotation>\n')
f.write(' <item key="Title" value="%s" sticky="true"/>\n' %title)
f.write(' <item key="AidaPath" value="%s" sticky="true"/>\n' %(path))
f.write(' <item key="FullPath" value="/%s.aida%s" sticky="true"/>\n' %(filename.split('/')[-1], path))
f.write(' </annotation>\n')
def write_datapointset_footer(self, f):
f.write(' </dataPointSet>\n')
def write_datapointset(self, f):
self.write_datapointset_header(f)
for bin, bindata in enumerate(self.data):
xval = 0.5*(bindata['UpEdge'] + bindata['LowEdge'])
if bindata['UpEdge'] == bindata['LowEdge']:
xerr = 0.5
else:
xerr = 0.5*(bindata['UpEdge'] - bindata['LowEdge'])
yval = bindata['Content']
yerr = bindata['Error']
self.write_datapoint(f, xval, xerr, yval, yerr)
self.write_datapointset_footer(f)
def remove_gaps(self):
# only look at histograms which are present in the reference file:
try:
refhist = refdata.histos['/REF%s' % self.description['AidaPath']]
except:
return
# check for differences in the binning and remove superfluous MC bins:
if len(refhist.data) != len(self.data):
print self.description['AidaPath']
newdata = []
for i in range(len(self.data)):
if self.data[i]['LowEdge'] == refhist.data[i]['LowEdge'] and \
self.data[i]['UpEdge'] == refhist.data[i]['UpEdge']:
newdata.append(self.data[i])
else:
print 'Deleted bin %d' %i
refhist.data.insert(i, self.data[i])
self.data = newdata
from optparse import OptionParser
parser = OptionParser(usage="%prog datafile MCfile outputfile")
opts, args = parser.parse_args()
if len(args) != 3:
sys.stderr.write("Must specity a reference, a MC, and an output file\n")
sys.exit(1)
# Convert the aida input files to flat files we can parse:
tempdir=tempfile.mkdtemp('.gap_removal')
filename = args[0].replace(".aida", "")
os.system("%s/aida2flat %s.aida > %s/%s.dat" %(os.path.dirname(os.path.realpath(sys.argv[0])), filename, tempdir, os.path.basename(filename)))
refdata = Inputdata("%s/%s" %(tempdir, os.path.basename(filename)))
filename = args[1].replace(".aida", "")
os.system("%s/aida2flat %s.aida > %s/%s.dat" %(os.path.dirname(os.path.realpath(sys.argv[0])), filename, tempdir, os.path.basename(filename)))
mcdata = Inputdata("%s/%s" %(tempdir, os.path.basename(filename)))
# Cleanup:
for i in os.listdir(tempdir):
os.unlink('%s/%s' %(tempdir, i))
os.rmdir(tempdir)
# Remove gap bins:
for i in mcdata.description['DrawOnly']:
mcdata.histos[i].remove_gaps()
# Write the new aida file with removed gap bins:
f = open(args[2], 'w')
f.write('<?xml version="1.0" encoding="ISO-8859-1" ?>\n')
f.write('<!DOCTYPE aida SYSTEM "http://aida.freehep.org/schemas/3.3/aida.dtd">\n')
f.write('<aida version="3.3">\n')
f.write(' <implementation version="1.1" package="FreeHEP"/>\n')
for i in mcdata.description['DrawOnly']:
mcdata.histos[i].write_datapointset(f)
f.write('</aida>\n')
f.close

File Metadata

Mime Type
text/x-python
Expires
Sat, Dec 21, 12:45 PM (1 d, 18 h)
Storage Engine
blob
Storage Format
Raw Data
Storage Handle
4022818
Default Alt Text
gap_removal (5 KB)

Event Timeline