5 import matplotlib.pyplot
as plt
15 KLP = KL + [
'rec.vtx.elastic.fuzzyk.png_idx']
19 df = tables[
'rec.vtx.elastic']
20 ret = (df[
'vtx.x'] > -180) & \
21 (df[
'vtx.x'] < 180) & \
22 (df[
'vtx.y'] > -180) & \
23 (df[
'vtx.y'] < 180) & \
27 kNDRockFilter =
Cut(kNDRockFilter)
31 dfvtx = tables[
'rec.vtx.elastic'][[
'vtx.x',
'vtx.y',
'vtx.z']]
32 dfpng = tables[
'rec.vtx.elastic.fuzzyk.png'][[
'start.x',
'start.y',
'start.z',
'rec.vtx.elastic.fuzzyk.png_idx']]
34 xdiff = dfvtx[
'vtx.x'] - dfpng[
'start.x']
35 ydiff = dfvtx[
'vtx.y'] - dfpng[
'start.y']
36 zdiff = dfvtx[
'vtx.z'] - dfpng[
'start.z']
37 length = np.sqrt(xdiff*xdiff + ydiff*ydiff + zdiff*zdiff)
39 idx = dfpng[
'rec.vtx.elastic.fuzzyk.png_idx']
40 ret = pd.concat([length, idx], axis=1)
41 ret = ret.reset_index().set_index(KLP)
43 return (ret[
'disp'] > 20).where(ret[
'disp'] > 20)
47 kProngNhitsCut =
Cut(
lambda tables: tables[
'rec.vtx.elastic.fuzzyk.png'][
'nhit'] < 6)
48 kPhotonID2018Cut =
Cut(
lambda tables: tables[
'rec.vtx.elastic.fuzzyk.png.cvnpart'][
'photonid'] < 0.6)
49 kIsNeutron =
Cut(
lambda tables:
abs(tables[
'rec.vtx.elastic.fuzzyk.png.truth'][
'motherpdg']) == 2112)
52 dfpng = tables[
'rec.vtx.elastic.fuzzyk.png.cvnpart'][
'muonid']
53 dflen = tables[
'rec.vtx.elastic.fuzzyk.png'][
'len']
55 dflen = dflen.reset_index().set_index(KL)
56 dfpng = dfpng.reset_index().set_index(KL)
57 maxid = dfpng[
'muonid'].groupby(level=KL).agg(np.max)
58 ret = pd.concat([dfpng[
'muonid'] - maxid, dfpng[
'rec.vtx.elastic.fuzzyk.png_idx'], dflen[
'len']], axis=1)
59 ret = ret.reset_index().set_index(KLP)
60 return ((ret[
'muonid'] < 0) & (ret[
'len'] <= 500))
63 kIgnoreMuon2018 =
Cut(kIgnoreMuon2018)
67 kProngCuts = kIgnoreMuon2018 & kIsNeutron & kProngNhitsCut & kPhotonID2018Cut
70 kProtonID =
Var(
lambda tables:
abs(tables[
'rec.vtx.elastic.fuzzyk.png.cvnpart'][
'protonid']))
73 filelist =
'/pnfs/nova/persistent/users/karlwarb/HDF5-Training-19-02-26/ND-ProngCVN-RHC/*.h5' 74 tablesProng = pa.loader(filelist,limit=limit, index=KLP)
75 tablesSlc = pa.loader(filelist, limit=limit)
77 a = pa.associate([tablesProng, tablesSlc])
80 specDummySlc = pa.spectrum(tablesSlc, kNDRockFilter & kNumuContainND, kProngDispl)
82 specProng = pa.spectrum(tablesProng, kProngCuts, kProtonID)
87 prongdispldf = specDummySlc.df()
88 vardf = specProng.df()
89 print(prongdispldf.head())
94 finaldf = pd.concat([prongdispldf, vardf], axis=1, join=
'inner')[vardf.name]
97 finalspec = pa.filledSpectrum(finaldf, specProng.POT())
100 n, bins = finalspec.histogram(bins=25, range=(0,1))
102 print(
'Selected ' +
str(n.sum()) +
' events from ' +
str(finalspec.POT()) +
' POT.')
104 plt.hist(bins[:-1], bins=bins, weights=n, histtype=
'step', color=
'blue', label=
'ND RHC')
105 plt.xlabel(
'Proton ID 2018')
108 plt.legend(loc=
'best')