Skip to content
Snippets Groups Projects
Commit f8ce75fb authored by Johannes Schumann's avatar Johannes Schumann
Browse files

Strip and update header

parent 90679798
No related branches found
No related tags found
1 merge request!6Detector geometries
Pipeline #19060 passed
......@@ -21,6 +21,7 @@ class DetectorVolume(ABC):
"""
def __init__(self):
self._volume = -1.0
self._coord_origin = (0., 0., 0.)
@abstractmethod
def random_pos(self):
......@@ -34,16 +35,39 @@ class DetectorVolume(ABC):
"""
pass
@abstractmethod
def header_entry(self):
"""
Returns the header information for the detector volume geometry
Returns
-------
tuple (header_key, header_information)
"""
pass
@property
def volume(self):
"""
Detector volume
Returns
-------
float [m^3]
The detector volume
"""
return self._volume
@property
def coord_origin(self):
"""
Coordinate origin
Returns
-------
tuple [m] (x, y, z)
"""
return self._coord_origin
class CanVolume(DetectorVolume):
"""
......@@ -76,6 +100,11 @@ class CanVolume(DetectorVolume):
pos_z = np.random.uniform(self._zmin, self._zmax)
return (pos_x, pos_y, pos_z)
def header_entry(self):
key = "can"
value = "{} {} {}".format(self._zmin, self._zmax, self._radius)
return key, value
class SphericalVolume(DetectorVolume):
"""
......@@ -89,10 +118,10 @@ class SphericalVolume(DetectorVolume):
Coordinate center of the sphere
(x, y, z)
"""
def __init__(self, radius, center=(0, 0, 0)):
def __init__(self, radius, coord_origin=(0, 0, 0)):
super().__init__()
self._radius = radius
self._center = center
self._coord_origin = coord_origin
self._volume = self._calc_volume()
def _calc_volume(self):
......@@ -106,4 +135,11 @@ class SphericalVolume(DetectorVolume):
pos_y = r * np.sin(phi) * np.sqrt(1 - np.power(cosTheta, 2))
pos_z = r * cosTheta
pos = (pos_x, pos_y, pos_z)
return tuple(np.add(self._center, pos))
return tuple(np.add(self._coord_origin, pos))
def header_entry(self):
key = "sphere"
value = "radius: {} center_x: {} center_y: {} center_z: {}".format(
self._radius, self._coord_origin[0], self._coord_origin[1],
self._coord_origin[2])
return key, value
......@@ -92,29 +92,16 @@ ROOTTUPLE_KEY = "RootTuple"
EMPTY_KM3NET_HEADER_DICT = {
"start_run": "0",
"XSecFile": "",
"drawing": "volume",
"detector": "",
"depth": "2475.0",
"muon_desc_file": "",
"target": "",
"cut_primary": "0 0 0 0",
"cut_seamuon": "0 0 0 0",
"cut_in": "0 0 0 0",
"cut_nu": "0 0 0 0",
"spectrum": "0",
"can": "0 0 0",
"flux": "0 0 0",
"fixedcan": "0 0 0 0 0",
"genvol": "0 0 0 0 0",
"coord_origin": "0 0 0",
"genhencut": "0 0",
"norma": "0 0",
"livetime": "0 0",
"seabottom": "0",
"DAQ": "0",
"tgen": "0",
"primary": "0",
"simul": ""
}
......@@ -441,8 +428,8 @@ def write_detector_file(gibuu_output,
media = read_default_media_compositions()
density = media["SeaWater"]["density"]
symbol = mendeleev.element(gibuu_output.Z).symbol
target = media["SeaWater"]["elements"][symbol]
element = mendeleev.element(gibuu_output.Z)
target = media["SeaWater"]["elements"][element.symbol]
target_density = 1e3 * density * target[1]
targets_per_volume = target_density * (1e3 * constants.Avogadro /
target[0].atomic_weight)
......@@ -452,16 +439,18 @@ def write_detector_file(gibuu_output,
head = ROOT.Head()
header_dct = EMPTY_KM3NET_HEADER_DICT.copy()
timestamp = datetime.now()
header_dct["simul"] = "KM3BUU {} {}".format(
version, timestamp.strftime("%Y%m%d %H%M%S"))
# header_dct["can"] = "{:.1f} {:.1f} {:.1f}".format(*can)
header_dct["tgen"] = "{:.1f}".format(livetime)
header_dct["target"] = element.name
key, value = geometry.header_entry()
header_dct[key] = value
header_dct["coord_origin"] = "{} {} {}".format(*geometry.coord_origin)
header_dct["flux"] = "{:d} 0 0".format(nu_type)
# header_dct["genvol"] = "0 {:.1f} {:.1f} {:.1f} {:d}".format(
# can[1], can[2], can_volume, gibuu_output.generated_events)
header_dct["cut_nu"] = "{:.2f} {:.2f} -1 1".format(gibuu_output.energy_min,
gibuu_output.energy_max)
header_dct["tgen"] = "{:.1f}".format(livetime)
header_dct["norma"] = "0 {}".format(gibuu_output.generated_events)
timestamp = datetime.now()
header_dct["simul"] = "KM3BUU {} {}".format(
version, timestamp.strftime("%Y%m%d %H%M%S"))
for k, v in header_dct.items():
head.set_line(k, v)
......
......@@ -20,8 +20,6 @@ class TestGeneralGeometry(unittest.TestCase):
with self.assertRaises(TypeError) as ctx:
d = DetectorVolume()
self.assertTrue('method random_pos' in str(ctx.exception))
class TestSphere(unittest.TestCase):
def setUp(self):
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment