1 | # $Id: extractShelve.py,v 1.3 2003/07/04 13:26:31 dressel Exp $ |
---|
2 | # ------------------------------------------------------------------- |
---|
3 | # GEANT4 tag $Name: $ |
---|
4 | # ------------------------------------------------------------------- |
---|
5 | # |
---|
6 | import string |
---|
7 | import math |
---|
8 | import os |
---|
9 | import shelve |
---|
10 | import string |
---|
11 | import detector |
---|
12 | |
---|
13 | |
---|
14 | |
---|
15 | # some functions for calculating the FOM |
---|
16 | |
---|
17 | def getR2(measure): |
---|
18 | "Error squared of a measure." |
---|
19 | v = measure.variance |
---|
20 | m2 = math.pow(measure.mean, 2) |
---|
21 | n = measure.entries |
---|
22 | |
---|
23 | r2 = 1e+10 |
---|
24 | if n>0 and m2 > 0: |
---|
25 | r2 = v/m2/n |
---|
26 | |
---|
27 | return r2 |
---|
28 | |
---|
29 | def getFOMbyName(shelveName, tallyName, bin): |
---|
30 | "FOM by shelve name, tally name and bin of energy region" |
---|
31 | she = shelve.open (shelveName, "r") |
---|
32 | measure = she[tallyName].measures[bin] |
---|
33 | R2= getR2(measure) |
---|
34 | T = she["runTime"] |
---|
35 | return 1./(T*R2) |
---|
36 | |
---|
37 | def getFOM(aShelve, tallyName, bin): |
---|
38 | "FOM given a shelve, tally name and energy bin." |
---|
39 | measure = aShelve[tallyName].measures[bin] |
---|
40 | R2= getR2(measure) |
---|
41 | T = aShelve["runTime"] |
---|
42 | return 1./(T*R2) |
---|
43 | |
---|
44 | |
---|
45 | |
---|
46 | |
---|
47 | |
---|
48 | |
---|
49 | |
---|
50 | |
---|
51 | |
---|
52 | |
---|
53 | # functions to retrive information from the shelve |
---|
54 | |
---|
55 | def infoShelve(file): |
---|
56 | "Print information stored in a shelve." |
---|
57 | she = shelve.open(file, "r") |
---|
58 | print file |
---|
59 | for k in she.keys(): |
---|
60 | value = she[k] |
---|
61 | if value.__class__.__module__ == "__builtin__" : |
---|
62 | print " ", k, ":", value |
---|
63 | |
---|
64 | def lsShelve(path): |
---|
65 | "List shelve information found in and under the given path. " |
---|
66 | if os.path.isfile(path): |
---|
67 | if string.find(path,".shelve") > -1: |
---|
68 | infoShelve(path) |
---|
69 | else: |
---|
70 | if os.path.isdir(path): |
---|
71 | files = os.listdir(path) |
---|
72 | for file in files: |
---|
73 | f=path + "/" + file |
---|
74 | lsShelve(f) |
---|
75 | |
---|
76 | |
---|
77 | |
---|
78 | def getFluxInRegion(she, dist="00", detType="ring", region=1): |
---|
79 | "Get the flux in a detector." |
---|
80 | tally = she["detector_" + dist + "Tally"] |
---|
81 | m = tally.measures[region] |
---|
82 | rawFlux = m.sum |
---|
83 | mean = m.mean |
---|
84 | var = m.getVariance() |
---|
85 | n = m.entries |
---|
86 | errRel = -1. |
---|
87 | if n>0 and var >0 and mean > 0: |
---|
88 | errRel = math.sqrt(var/n)/mean |
---|
89 | nPeakNeutrons = she["generatorTally"].measures[1].sum |
---|
90 | scale = detector.detScale(nPeakNeutrons, she["energy"], 0) |
---|
91 | return rawFlux * scale / detector.DetectorVolume[detType][dist], errRel |
---|
92 | |
---|
93 | |
---|
94 | |
---|
95 | def getFluxes(she, detType = "ring"): |
---|
96 | "Get all fluxes and FOM in a shelve." |
---|
97 | fluxName = "flux_" + she["energy"]+"_"+\ |
---|
98 | she["shieldWidth"] |
---|
99 | keys = she.keys() |
---|
100 | fluxes = {} |
---|
101 | for k in keys: |
---|
102 | if string.find(k,"detector_") > -1: |
---|
103 | sp = string.split(k,"detector_") |
---|
104 | sp = string.split(sp[1],"Tally") |
---|
105 | dist = sp[0] |
---|
106 | fluxNameD = fluxName + "_" + dist |
---|
107 | for i in range(1,3): |
---|
108 | fom = getFOM(she, k, i) |
---|
109 | f = fluxNameD + "_" + "%(i)d" % vars() |
---|
110 | flux, errRel = getFluxInRegion(she, dist, detType, i) |
---|
111 | errRel*=100 |
---|
112 | sFlux = '%1.2E' % (flux) |
---|
113 | sErrRel = '%1.0f' % (errRel) |
---|
114 | if errRel < 10.0: |
---|
115 | sErrRel = '%1.1f' % (errRel) |
---|
116 | print f, sFlux, ", errRel:", sErrRel, " FOM:", fom |
---|
117 | # fluxes[f] = sFlux |
---|
118 | return fluxes |
---|