source: trunk/environments/g4py/site-modules/utils/MCScore/mcscore.py @ 1337

Last change on this file since 1337 was 1337, checked in by garnier, 14 years ago

tag geant4.9.4 beta 1 + modifs locales

File size: 3.5 KB
Line 
1"""
2Python module
3
4This module provides classes and functions for scoring reactions
5
6  [C] MCVertex:
7  [C] MCParticle:
8  [f] read_next_vertex(stream):
9
10                                              Q, 2006
11"""
12import string
13from Geant4.hepunit import *
14
15# ==================================================================
16# public symbols
17# ==================================================================
18__all__ = [ 'MCParticle', 'MCVertex', 'read_next_vertex' ]
19
20
21# ==================================================================
22#   class definition
23# ==================================================================
24
25# ------------------------------------------------------------------
26# MCParticle
27# ------------------------------------------------------------------
28class MCParticle:
29  "MC particle"
30  def __init__(self, aname, aZ, aA, akE, apx, apy, apz):
31    self.name = aname
32    self.Z = aZ
33    self.A = aA
34    self.kineticE = akE
35    self.px = apx
36    self.py = apy
37    self.pz = apz
38
39  def printout(self):
40    print "--- particle: %s, Z=%2d, A=%2d, kE=%g" % \
41          (self.name, self.Z, self.A, self.kineticE/MeV)
42
43# ------------------------------------------------------------------
44# MCVertex
45# ------------------------------------------------------------------
46class MCVertex :
47  "MC vertex"
48  def __init__(self, ax, ay, az):
49    self.x = ax
50    self.y = ay
51    self.z = az
52    self.nparticle = 0
53    self.particle_list = []
54
55  def append_particle(self, aparticle):
56    self.particle_list.append(aparticle)
57    self.nparticle= self.nparticle+1
58
59  def printout(self):
60    print "@@@ vertex: x=(%g,%g,%g) Nsec=%3d" % \
61          (self.x/cm, self.y/cm, self.z/cm, self.nparticle)
62    for p in self.particle_list:
63      p.printout()
64 
65  def dump_vertex(self, stream):
66    aline = "%g %g %g %d\n" % \
67            (self.x/m, self.y/m, self.z/m, self.nparticle)
68    stream.write(aline)
69    for p in self.particle_list:
70      aline = " %s %d %d %g %g %g %g\n" % \
71              (p.name, p.Z, p.A, p.kineticE/MeV, p.px/MeV, p.py/MeV, p.pz/MeV)
72      stream.write(aline)
73
74  def __del__(self):
75    np = len(self.particle_list)
76    del self.particle_list[0:np]
77   
78
79# ==================================================================
80#   I/O interface
81# ==================================================================
82def read_next_vertex(stream):
83  "read next vertex from a file stream"
84  line= stream.readline()
85  if line == "":  # EOF
86    return 0
87
88  # reading vertex
89  data = line.split()
90  x = string.atof(data[0]) * m
91  y = string.atof(data[1]) * m
92  z = string.atof(data[2]) * m
93  nsec = string.atoi(data[3])
94
95  vertex = MCVertex(x,y,z)
96
97  # reading particles
98  for p in range(0, nsec):
99    data = stream.readline().split()
100    pname = data[0]
101    Z = string.atoi(data[1])
102    A = string.atoi(data[2])   
103    kE = string.atof(data[3]) * MeV
104    px = string.atof(data[4]) * MeV
105    py = string.atof(data[5]) * MeV
106    pz = string.atof(data[6]) * MeV
107
108    particle = MCParticle(pname, Z, A, kE, px, py, pz)
109    vertex.append_particle(particle)
110
111  return vertex
112
113
114# ==================================================================
115# test
116# ==================================================================
117def test():
118  f = open("reaction.dat")
119  f.seek(0)
120
121  while(1):
122    vertex = read_next_vertex(f)
123    if vertex == 0:
124      break
125    vertex.printout()
126    del vertex
127  f.close()
128  print ">>> EOF"
129
130
131# ==================================================================
132# main
133# ==================================================================
134if __name__ == "__main__":
135  test()
136 
Note: See TracBrowser for help on using the repository browser.