Belle II Software  release-08-01-10
GeoARICHBtestCreator.cc
1 /**************************************************************************
2  * basf2 (Belle II Analysis Software Framework) *
3  * Author: The Belle II Collaboration *
4  * *
5  * See git log for contributors and copyright holders. *
6  * This file is licensed under LGPL-3.0, see LICENSE.md. *
7  **************************************************************************/
8 #include <sstream>
9 #include <string.h>
10 #include <boost/format.hpp>
11 #include <boost/foreach.hpp>
12 #include <boost/algorithm/string.hpp>
13 
14 //#include <geant4/G4LogicalVolume.hh>
15 
16 // Geant4
17 #include <G4LogicalVolume.hh>
18 #include <G4PVPlacement.hh>
19 #include <G4LogicalSkinSurface.hh>
20 #include <G4OpticalSurface.hh>
21 // Geant4 Shapes
22 #include <G4Box.hh>
23 #include <G4SubtractionSolid.hh>
24 #include <G4Material.hh>
25 
26 #include <Python.h>
27 
28 #include <arich/geometry/GeoARICHBtestCreator.h>
29 #include <arich/geometry/ARICHGeometryPar.h>
30 #include <arich/geometry/ARICHBtestGeometryPar.h>
31 
32 #include <geometry/Materials.h>
33 #include <geometry/CreatorFactory.h>
34 #include <geometry/utilities.h>
35 #include <framework/gearbox/GearDir.h>
36 #include <framework/gearbox/Unit.h>
37 #include <framework/logging/Logger.h>
38 // Framework - DataStore
39 #include <framework/datastore/StoreObjPtr.h>
40 #include <framework/dataobjects/EventMetaData.h>
41 
42 #include <arich/simulation/SensitiveDetector.h>
43 #include <arich/simulation/SensitiveAero.h>
44 
45 using namespace std;
46 using namespace boost;
47 
48 namespace Belle2 {
54  using namespace geometry;
55 
56  namespace arich {
57 
58  //-----------------------------------------------------------------
59  // Register the Creator
60  //-----------------------------------------------------------------
61 
62  geometry::CreatorFactory<GeoARICHBtestCreator> GeoARICHBtestFactory("ARICHBtestCreator");
63 
64  //-----------------------------------------------------------------
65  // Implementation
66  //-----------------------------------------------------------------
67 
68  GeoARICHBtestCreator::GeoARICHBtestCreator():
69  m_sensitive(new SensitiveDetector),
70  m_sensitiveAero(new SensitiveAero),
71  m_runno(0),
72  m_neve(0),
73  m_rotation(0),
74  m_rx(0),
75  m_ry(0),
76  m_aerosupport(0),
77  m_aerogeldx(0),
78  m_framedx(0),
79  m_rotation1(0),
80  m_configuration(0)
81  {};
82 
83  GeoARICHBtestCreator::~GeoARICHBtestCreator()
84  {
85 
86  }
87 
88 
89  void GeoARICHBtestCreator::create(const GearDir& content, G4LogicalVolume& topVolume, GeometryTypes)
90  {
91 
92  B2INFO("GeoARICHBtestCreator::create");
93  StoreObjPtr<EventMetaData> eventMetaDataPtr;
94 
95  int run = 68;
96  PyObject* m = PyImport_AddModule(strdup("__main__"));
97  if (m) {
98  PyObject* v = PyObject_GetAttrString(m, strdup("runno"));
99  if (v) {
100  run = PyLong_AsLong(v);
101  Py_DECREF(v);
102  }
103  B2INFO("GeoARICHBtestCreator::create runno = " << run);
104  }
105 
106  B2INFO("eventMetaDataPtr run:" << run);
107  // eventMetaDataPtr->setEndOfData();
108 
109 
110 
111 
112  string Type = content.getString("@type", "");
113 
114  char nodestr[100];
115  sprintf(nodestr, "run[runno=%d]", run);
116  if (Type == "beamtest") {
117  BOOST_FOREACH(const GearDir & runparam, content.getNodes(nodestr)) {
118  m_runno = runparam.getInt("runno", -1);
119  m_author = runparam.getString("author", "");
120  m_neve = runparam.getInt("neve", -1);
121  m_runtype = runparam.getString("calibration", "pion");
122  m_hapdID = runparam.getString("setup1", "unknown");
123  m_aerogelID = runparam.getString("aerogel1", "unknown");
124  m_mirrorID = runparam.getString("mirror", "unknown");
125  m_rotation = runparam.getDouble("rotation", 0);
126  m_rx = runparam.getDouble("positionx", 0);
127  m_ry = runparam.getDouble("positiony", 0);
128  m_mytype = runparam.getString("type1", "unknown");
129  m_daqqa = runparam.getString("daqqa1", "unknown");
130  m_comment = runparam.getString("comment1", "unknown");
131  m_datum = runparam.getString("datum", "unknown");
132 
133 
134  B2INFO("runno : " << m_runno);
135  B2INFO("author : " << m_author);
136  B2INFO("neve : " << m_neve);
137  B2INFO("runtype : " << m_runtype);
138  B2INFO("hapdID : " << m_hapdID);
139  B2INFO("aerogelID: " << m_aerogelID);
140  B2INFO("mirrorID : " << m_mirrorID);
141  B2INFO("rotation : " << m_rotation);
142  B2INFO("rx : " << m_rx);
143  B2INFO("ry : " << m_ry);
144  B2INFO("runtype : " << m_mytype);
145  B2INFO("daqqa : " << m_daqqa);
146  B2INFO("comment : " << m_comment);
147  B2INFO("datum : " << m_datum);
148 
149 
150  }
151  string aerogelname;
152  sprintf(nodestr, "setup/aerogel/row[@id=\"%s\"]", m_aerogelID.c_str());
153 
154  GearDir runparam(content, nodestr);
155  B2INFO("id : " << runparam.getString("@id", ""));
156  BOOST_FOREACH(const GearDir & aeroparam, runparam.getNodes("aerogel")) {
157  aerogelname = aeroparam.getString(".", "");
158  string stype = aeroparam.getString("@type", "");
159  B2INFO(stype << " aerogelname : " << aerogelname);
160  sprintf(nodestr, "setup/aerogelinfo/row[@id=\"%s\"]", aerogelname.c_str());
161  GearDir infoparam(content, nodestr);
162 
163  double agelrefind = infoparam.getDouble("refind", 1);
164  double ageltrlen = infoparam.getLength("trlen", 0);
165  double agelthickness = infoparam.getLength("thickness", 0);
166  if (stype != string("left")) {
167  m_ageltrlen.push_back(ageltrlen);
168  m_agelrefind.push_back(agelrefind);
169  m_agelthickness.push_back(agelthickness);
170  }
171  B2INFO("refind : " << agelrefind);
172  B2INFO("trlen : " << ageltrlen / Unit::mm);
173  B2INFO("thickness : " << agelthickness / Unit::mm);
174 
175  }
176  int size = m_hapdID.size();
177 
178  m_aerosupport = 0;
179  if (size > 0) {
180  char agelsupport = m_hapdID.at(size - 1);
181  if (agelsupport == 'a') m_aerosupport = 1;
182  if (agelsupport == 'b') m_aerosupport = 2;
183  }
184 
185  if (m_aerosupport) size--;
186  sprintf(nodestr, "setup/hapd/row[@id=\"%s\"]", m_hapdID.substr(0, size).c_str());
187  B2INFO("nodestr : " << nodestr);
188  B2INFO("aerogelsupport : " << m_aerosupport);
189  GearDir hapdparam(content, nodestr);
190  //BOOST_FOREACH(const GearDir & runparam, content.getNodes(nodestr)) {
191  m_aerogeldx = hapdparam.getLength("aerogeldx", 0);
192  m_framedx = hapdparam.getLength("framedx", 0) * CLHEP::mm / Unit::mm ;
193  m_rotation1 = hapdparam.getDouble("rotation", 0);
194  m_configuration = hapdparam.getInt("setup", 0);
195  m_comment1 = hapdparam.getString("comment", "");
196 
197 
198  B2INFO("aerogeldx : " << m_aerogeldx);
199  B2INFO("framedx : " << m_framedx);
200  B2INFO("rotation : " << m_rotation1);
201  B2INFO("configuration : " << m_configuration);
202  B2INFO("comment : " << m_comment);
203  //}
204 
205  GearDir setup(content, "setup");
206 
207  createBtestGeometry(setup, topVolume);
208  }
209  }
210 
211  double GeoARICHBtestCreator::getAvgRINDEX(G4Material* material)
212  {
213  G4MaterialPropertiesTable* mTable = material->GetMaterialPropertiesTable();
214  if (!mTable) return 0;
215  G4MaterialPropertyVector* mVector = mTable->GetProperty("RINDEX");
216  if (!mVector) return 0;
217  G4bool b;
218  return mVector->GetValue(2 * Unit::eV / Unit::MeV, b);
219  }
220 
221  G4LogicalVolume* GeoARICHBtestCreator::buildModule(GearDir Module)
222  {
223 
224  // get detector module parameters
225 
226  // get module materials
227  string wallMat = Module.getString("wallMaterial");
228  string winMat = Module.getString("windowMaterial");
229  string botMat = Module.getString("Bottom/material");
230  G4Material* wallMaterial = Materials::get(wallMat);
231  G4Material* windowMaterial = Materials::get(winMat);
232  G4Material* bottomMaterial = Materials::get(botMat);
233  G4Material* boxFill = Materials::get("ARICH_Vacuum");
234 
235  // check that module window material has specified refractive index
236  double wref = getAvgRINDEX(windowMaterial);
237  if (!wref) B2WARNING("Material '" << winMat <<
238  "', required for ARICH photon detector window as no specified refractive index. Continuing, but no photons in ARICH will be detected.");
239  ARICHGeometryPar* m_arichgp = ARICHGeometryPar::Instance();
240  m_arichgp->setWindowRefIndex(wref);
241  // get module dimensions
242  double modXsize = Module.getLength("moduleXSize") / Unit::mm;
243  double modZsize = Module.getLength("moduleZSize") / Unit::mm;
244  double wallThick = Module.getLength("moduleWallThickness") / Unit::mm;
245  double winThick = Module.getLength("windowThickness") / Unit::mm ;
246  double sensXsize = m_arichgp->getSensitiveSurfaceSize() / Unit::mm;
247  double botThick = Module.getLength("Bottom/thickness") / Unit::mm;
248 
249  // some trivial checks of overlaps
250  if (sensXsize > modXsize - 2 * wallThick)
251  B2FATAL("ARICH photon detector module: Sensitive surface is too big. Doesn't fit into module box.");
252  if (winThick + botThick > modZsize)
253  B2FATAL("ARICH photon detector module: window + bottom thickness larger than module thickness.");
254 
255  // module master volume
256  G4Box* moduleBox = new G4Box("Box", modXsize / 2., modXsize / 2., modZsize / 2.);
257  G4LogicalVolume* lmoduleBox = new G4LogicalVolume(moduleBox, boxFill, "moduleBox");
258 
259  // build and place module wall
260  G4Box* tempBox = new G4Box("tempBox", modXsize / 2. - wallThick, modXsize / 2. - wallThick,
261  modZsize / 2. + 0.1); // Dont't care about "+0.1", needs to be there.
262  G4SubtractionSolid* moduleWall = new G4SubtractionSolid("Box-tempBox", moduleBox, tempBox);
263  G4LogicalVolume* lmoduleWall = new G4LogicalVolume(moduleWall, wallMaterial, "moduleWall");
264  setColor(*lmoduleWall, "rgb(1.0,0.0,0.0,1.0)");
265  new G4PVPlacement(G4Transform3D(), lmoduleWall, "moduleWall", lmoduleBox, false, 1);
266 
267  // build module window
268  G4Box* winBox = new G4Box("winBox", modXsize / 2. - wallThick, modXsize / 2. - wallThick, winThick / 2.);
269  G4LogicalVolume* lmoduleWin = new G4LogicalVolume(winBox, windowMaterial, "moduleWindow");
270  setColor(*lmoduleWin, "rgb(0.7,0.7,0.7,1.0)");
271  G4Transform3D transform = G4Translate3D(0., 0., (-modZsize + winThick) / 2.);
272  new G4PVPlacement(transform, lmoduleWin, "moduleWindow", lmoduleBox, false, 1);
273 
274  // build module bottom
275  G4Box* botBox = new G4Box("botBox", modXsize / 2. - wallThick, modXsize / 2. - wallThick, botThick / 2.);
276  G4LogicalVolume* lmoduleBot = new G4LogicalVolume(botBox, bottomMaterial, "moduleBottom");
277  // if (isBeamBkgStudy) lmoduleBot->SetSensitiveDetector(new BkgSensitiveDetector("ARICH", 1));
278  setColor(*lmoduleBot, "rgb(0.0,1.0,0.0,1.0)");
279  G4Transform3D transform1 = G4Translate3D(0., 0., (modZsize - botThick) / 2.);
280  // add surface optical properties if specified
281  Materials& materials = Materials::getInstance();
282  GearDir bottomParam(Module, "Bottom/Surface");
283  if (bottomParam) {
284  G4OpticalSurface* optSurf = materials.createOpticalSurface(bottomParam);
285  new G4LogicalSkinSurface("bottomSurface", lmoduleBot, optSurf);
286  } else B2INFO("ARICH: No optical properties are specified for detector module bottom surface.");
287  new G4PVPlacement(transform1, lmoduleBot, "moduleBottom", lmoduleBox, false, 1);
288 
289  // build sensitive surface
290  G4Box* sensBox = new G4Box("sensBox", sensXsize / 2., sensXsize / 2., 0.1 * Unit::mm);
291  G4LogicalVolume* lmoduleSens = new G4LogicalVolume(sensBox, boxFill, "moduleSensitive");
292  lmoduleSens->SetSensitiveDetector(m_sensitive);
293  setColor(*lmoduleSens, "rgb(0.5,0.5,0.5,1.0)");
294  G4Transform3D transform2 = G4Translate3D(0., 0., (-modZsize + 0.1) / 2. + winThick);
295  new G4PVPlacement(transform2, lmoduleSens, "moduleSensitive", lmoduleBox, false, 1);
296 
297  // module is build, return module logical volume
298  return lmoduleBox;
299  }
300 
301 
302  G4Material* GeoARICHBtestCreator::createAerogel(const char* aeroname, double RefractiveIndex, double AerogelTransmissionLength)
303  {
304 
305 
306  G4double density = (RefractiveIndex - 1) / 0.21 * CLHEP::g / CLHEP::cm3;
307  B2INFO("Creating ARICH " << aeroname << " n=" << RefractiveIndex << " density=" << density / CLHEP::g * CLHEP::cm3 << " g/cm3");
308  Materials& materials = Materials::getInstance();
309  G4Material* _aerogel = new G4Material(aeroname, density, 4);
310  _aerogel->AddElement(materials.getElement("O"), 0.665);
311  _aerogel->AddElement(materials.getElement("H"), 0.042);
312  _aerogel->AddElement(materials.getElement("Si"), 0.292);
313  _aerogel->AddElement(materials.getElement("C"), 0.001);
314 
315 
316  const G4double AerogelAbsorbtionLength = 1000 * Unit::mm;
317 
318  const G4int NBins = 40;
319  G4double MomentumBins[NBins];
320 
321  G4double AerogelRindex[NBins];
322  G4double AerogelAbsorption[NBins];
323  G4double AerogelRayleigh[NBins];
324 
325  G4double MaxPhotonEnergy = 5 * CLHEP::eV;
326  G4double MinPhotonEnergy = 1.5 * CLHEP::eV;
327 
328  for (G4int i = 0; i < NBins; i++) {
329 
330  const G4double energy = float(i) / NBins * (MaxPhotonEnergy - MinPhotonEnergy) + MinPhotonEnergy;
331 
332  MomentumBins[i] = energy;
333  AerogelRindex[i] = RefractiveIndex;
334  AerogelAbsorption[i] = AerogelAbsorbtionLength;
335 
336  const G4double Lambda0 = 400 * 1e-6 * CLHEP::mm;
337  const G4double Lambda = 1240 * CLHEP::eV / energy * 1e-6 * CLHEP::mm;
338  G4double x = Lambda / Lambda0;
339  AerogelRayleigh[i] = AerogelTransmissionLength * x * x * x * x;
340  }
341 
342 
343  G4MaterialPropertiesTable* AeroProperty = new G4MaterialPropertiesTable();
344  AeroProperty->AddProperty("RINDEX", MomentumBins, AerogelRindex, NBins);
345  AeroProperty->AddProperty("ABSLENGTH", MomentumBins, AerogelAbsorption, NBins);
346  AeroProperty->AddProperty("RAYLEIGH", MomentumBins, AerogelRayleigh, NBins);
347 
348 
349  _aerogel->SetMaterialPropertiesTable(AeroProperty);
350 
351 
352  return _aerogel;
353  }
354 
355 
356  void GeoARICHBtestCreator::createBtestGeometry(const GearDir& content, G4LogicalVolume& topWorld)
357  {
358 
359  B2INFO("ARICH Btest geometry will be built.");
360  ARICHGeometryPar* m_arichgp = ARICHGeometryPar::Instance();
361 
362  ARICHBtestGeometryPar* m_arichbtgp = ARICHBtestGeometryPar::Instance();
363 
364  // experimental box
365 
366  GearDir boxParams(content, "ExperimentalBox");
367  double xBox = boxParams.getLength("xSize") * CLHEP::mm / Unit::mm;
368  double yBox = boxParams.getLength("ySize") * CLHEP::mm / Unit::mm;
369  double zBox = boxParams.getLength("zSize") * CLHEP::mm / Unit::mm;
370 
371  double xoffset = boxParams.getLength("beamcenter/x") * CLHEP::mm / Unit::mm;
372  double yoffset = boxParams.getLength("beamcenter/y") * CLHEP::mm / Unit::mm;
373  double zoffset = boxParams.getLength("beamcenter/z") * CLHEP::mm / Unit::mm - zBox / 2.;
374  G4ThreeVector roffset(xoffset, yoffset, zoffset);
375 
376  TVector3 sh(boxParams.getLength("beamcenter/x"), boxParams.getLength("beamcenter/y"),
377  boxParams.getLength("beamcenter/z") - boxParams.getLength("zSize") / 2.);
378  m_arichbtgp->setOffset(sh);
379 
380  string boxMat = boxParams.getString("material");
381  G4Material* boxMaterial = Materials::get(boxMat);
382  G4Box* expBox = new G4Box("ExperimentalBox", xBox / 2., yBox / 2., zBox / 2.);
383  G4LogicalVolume* topVolume = new G4LogicalVolume(expBox, boxMaterial, "ARICH.experimentalbox");
384  new G4PVPlacement(G4Transform3D(), topVolume, "ARICH.experimentalbox", &topWorld, false, 1);
385  setVisibility(*topVolume, false);
386 
387  TVector3 trackingshift(content.getLength("tracking/shift/x"),
388  content.getLength("tracking/shift/y"),
389  content.getLength("tracking/shift/z"));
390 
391  char mnodestr[256];
392  sprintf(mnodestr, "tracking/shift/run[@id=\"%d\"]", m_runno);
393  if (content.exists(mnodestr)) {
394  GearDir runtrackingshift(content, mnodestr);
395  trackingshift[0] = runtrackingshift.getLength("x");
396  trackingshift[1] = runtrackingshift.getLength("y");
397  trackingshift[2] = runtrackingshift.getLength("z");
398  }
399  m_arichbtgp->setTrackingShift(trackingshift);
400  ARICHTracking* m_mwpc = new ARICHTracking[4];
401  m_arichbtgp->setMwpc(m_mwpc);
402  BOOST_FOREACH(const GearDir & mwpc, content.getNodes("tracking/mwpc")) {
403  double x = mwpc.getLength("size/x") * CLHEP::mm / Unit::mm;
404  double y = mwpc.getLength("size/y") * CLHEP::mm / Unit::mm;
405  double z = mwpc.getLength("size/z") * CLHEP::mm / Unit::mm;
406 
407  double px = mwpc.getLength("position/x") * CLHEP::mm / Unit::mm;
408  double py = mwpc.getLength("position/y") * CLHEP::mm / Unit::mm;
409  double pz = mwpc.getLength("position/z") * CLHEP::mm / Unit::mm;
410 
411  G4Box* mwpcBox = new G4Box("MwpcBox", x / 2., y / 2., z / 2.);
412  G4LogicalVolume* mwpcVol = new G4LogicalVolume(mwpcBox, Materials::get(mwpc.getString("material")), "ARICH.mwpc");
413  new G4PVPlacement(G4Transform3D(G4RotationMatrix(), G4ThreeVector(px, py, pz) + roffset), mwpcVol, "ARICH.mwpc", topVolume, false,
414  1);
415  //setVisibility(*mwpc, true);
416 
417  int id = mwpc.getInt("@id", -1);
418  B2INFO("GeoARICHBtestCreator::" << LogVar("MWPC ID", id));
419  if (id < 4 && id >= 0) {
420  m_mwpc[id].tdc[0] = mwpc.getInt("tdc/y/up");
421  m_mwpc[id].tdc[1] = mwpc.getInt("tdc/y/down");
422  m_mwpc[id].tdc[2] = mwpc.getInt("tdc/x/left");
423  m_mwpc[id].tdc[3] = mwpc.getInt("tdc/x/right");
424  m_mwpc[id].atdc = mwpc.getInt("tdc/anode", 0);
425  m_mwpc[id].slp[0] = mwpc.getDouble("slope/x");
426  m_mwpc[id].slp[1] = mwpc.getDouble("slope/y");
427  m_mwpc[id].offset[0] = mwpc.getDouble("offset/x");
428  m_mwpc[id].offset[1] = mwpc.getDouble("offset/y");
429  m_mwpc[id].cutll[0] = mwpc.getInt("tdccut/y/min");
430  m_mwpc[id].cutll[1] = mwpc.getInt("tdccut/x/min");
431  m_mwpc[id].cutul[0] = mwpc.getInt("tdccut/y/max");
432  m_mwpc[id].cutul[1] = mwpc.getInt("tdccut/x/max");
433  m_mwpc[id].pos[0] = mwpc.getDouble("position/x");
434  m_mwpc[id].pos[1] = mwpc.getDouble("position/y");
435  m_mwpc[id].pos[2] = mwpc.getDouble("position/z");
436  // m_mwpc[id].Print();
437  }
438 
439  }
440  // physical position of the hapd channels
441 
442  istringstream mapstream;
443  double mx, my;
444  mapstream.str(content.getString("hapdmap"));
445  while (mapstream >> mx >> my) {
446  m_arichbtgp->AddHapdChannelPositionPair(mx, my);
447  }
448  mapstream.clear();
449 
450  // mapping of the electronic channels
451  int ipx, ipy;
452  mapstream.str(content.getString("hapdchmap"));
453  while (mapstream >> ipx >> ipy) {
454  m_arichbtgp->AddHapdElectronicMapPair(ipx, ipy);
455  }
456  // experimental frame consisting of detector plane, aerogel and mirrors
457 
458  GearDir frameParams(content, "Frame");
459  double xFrame = frameParams.getLength("xSize") * CLHEP::mm / Unit::mm;
460  double yFrame = frameParams.getLength("ySize") * CLHEP::mm / Unit::mm;
461  double zFrame = frameParams.getLength("zSize") * CLHEP::mm / Unit::mm;
462  string envMat = frameParams.getString("material");
463 
464  double px = frameParams.getLength("position/x") * CLHEP::mm / Unit::mm;
465  double py = frameParams.getLength("position/y") * CLHEP::mm / Unit::mm;
466  double pz = frameParams.getLength("position/z") * CLHEP::mm / Unit::mm;
467 
468  G4Material* envMaterial = Materials::get(envMat);
469 
470 
471  G4Box* envBox = new G4Box("FrameBox", xFrame / 2., yFrame / 2., zFrame / 2.);
472  G4LogicalVolume* lenvBox = new G4LogicalVolume(envBox, envMaterial, "ARICH.frame");
473  G4ThreeVector frameOrigin0(m_framedx + px, py, pz); // rotation point of the detector frame wrt beamcenter
474  G4ThreeVector frameOrigin = frameOrigin0 + roffset;
475  G4RotationMatrix frameRotation;
476  frameRotation.rotateY(-m_rotation1 * CLHEP::degree);
477  G4Transform3D frameTransformation = G4Transform3D(frameRotation, frameOrigin);
478 
479  new G4PVPlacement(frameTransformation, lenvBox, "ARICH.frame", topVolume, false, 1);
480  //setVisibility(*lenvBox, false);
481 
482  TVector3 rotationCenter = TVector3(frameOrigin0.x() * Unit::mm / CLHEP::mm, frameOrigin0.y() * Unit::mm / CLHEP::mm,
483  frameOrigin0.z() * Unit::mm / CLHEP::mm);
484  m_arichbtgp->setFrameRotation(m_rotation1 * CLHEP::degree);
485  m_arichbtgp->setRotationCenter(rotationCenter);
486 
487 
488  char nodestr[256];
489  B2INFO(content.getPath());
490  sprintf(nodestr, "PhotonDetector/setup[@id=\"%d\"]", m_configuration);
491  GearDir hapdcontent(content, nodestr);
492  B2INFO(hapdcontent.getPath());
493 
494 
495 
496  char mirrornodestr[256];
497  sprintf(mirrornodestr, "Mirrors/setup[@id=\"%s\"]", m_mirrorID.c_str());
498 
499  GearDir mirrorcontent(content, mirrornodestr);
500  B2INFO(mirrorcontent.getPath());
501 
502  // detectors
503  m_arichgp->Initialize(hapdcontent, mirrorcontent);
504 
505 
506  GearDir moduleParam(hapdcontent, "Detector/Module");
507  G4LogicalVolume* detModule = buildModule(moduleParam);
508 
509  double detZpos = hapdcontent.getLength("Detector/Plane/zPosition") * CLHEP::mm / Unit::mm;
510  double detThick = hapdcontent.getLength("Detector/Module/moduleZSize") * CLHEP::mm / Unit::mm;
511  int nModules = m_arichgp->getNMCopies();
512 
513  for (int i = 1; i <= nModules; i++) {
514  G4ThreeVector origin = m_arichgp->getOriginG4(i);
515  origin.setZ(detZpos + detThick / 2.);
516  double angle = m_arichgp->getModAngle(i);
517  G4RotationMatrix Ra;
518  Ra.rotateZ(angle);
519  G4Transform3D trans = G4Transform3D(Ra, origin);
520  new G4PVPlacement(G4Transform3D(Ra, origin), detModule, "detModule", lenvBox, false, i);
521  B2INFO(nodestr << "Module " << i << " is build ");
522  }
523  // mask hot channels
524  int npx = m_arichgp->getDetectorXPadNumber();
525  BOOST_FOREACH(const double & ch, hapdcontent.getArray("HotChannels")) {
526  int channelID = (int) ch;
527  int moduleID = (npx) ? channelID / (npx * npx) : 0;
528  channelID %= (npx * npx);
529  m_arichgp->setActive(moduleID, channelID, false);
530  B2INFO("HotChannel " << ch << " : Module " << moduleID << "channelID " << channelID << " disabled");
531  }
532  // mask dead channels
533  BOOST_FOREACH(const double & ch, hapdcontent.getArray("DeadChannels")) {
534  int channelID = (int) ch;
535  int moduleID = (npx) ? channelID / (npx * npx) : 0;
536  channelID %= (npx * npx);
537  m_arichgp->setActive(moduleID, channelID, false);
538  B2INFO("DeadChannel " << ch << " : Module " << moduleID << "channelID " << channelID << " disabled");
539  }
540  // place aerogel tiles
541  GearDir aerogelParam(content, "Aerogel");
542  double sizeX = aerogelParam.getLength("tileXSize") * CLHEP::mm / Unit::mm;
543  double sizeY = aerogelParam.getLength("tileYSize") * CLHEP::mm / Unit::mm;
544  double posX = aerogelParam.getLength("tileXPos") * CLHEP::mm / Unit::mm;
545  double posY = aerogelParam.getLength("tileYPos") * CLHEP::mm / Unit::mm;
546  double posZ = aerogelParam.getLength("tileZPos") * CLHEP::mm / Unit::mm;
547  double posZ0 = posZ;
548  double meanrefind = 0;
549  double meantrlen = 0;
550 
551  // get parameter from python script
552  PyObject* m = PyImport_AddModule(strdup("__main__"));
553  if (m) {
554  int averageagel = 0;
555  PyObject* v = PyObject_GetAttrString(m, strdup("averageagel"));
556  if (v) {
557  averageagel = PyLong_AsLong(v);
558  Py_DECREF(v);
559  }
560  B2INFO("Python averageagel = " << averageagel);
561  m_arichbtgp->setAverageAgel(averageagel > 0);
562  }
563 
564  for (unsigned int ilayer = 0; ilayer < m_agelthickness.size(); ilayer++) {
565  char aeroname[100];
566  sprintf(aeroname, "Aerogel%u", ilayer + 1);
567  G4Material* tileMaterial = createAerogel(aeroname, m_agelrefind[ilayer], m_ageltrlen[ilayer]);
568  double sizeZ = m_agelthickness[ilayer] * CLHEP::mm / Unit::mm;
569 
570  if (!m_arichbtgp->getAverageAgel()) {
571  m_arichgp->setAeroRefIndex(ilayer, m_agelrefind[ilayer]);
572  m_arichgp->setAerogelZPosition(ilayer, (posZ - zFrame / 2.) * Unit::mm / CLHEP::mm);
573  m_arichgp->setAerogelThickness(ilayer, sizeZ * Unit::mm / CLHEP::mm);
574  m_arichgp->setAeroTransLength(ilayer, m_ageltrlen[ilayer]);
575  }
576 
577  meantrlen += sizeZ / m_ageltrlen[ilayer];
578  meanrefind += m_agelrefind[ilayer];
579  G4Box* tileBox = new G4Box("tileBox", sizeX / 2., sizeY / 2., sizeZ / 2.);
580  G4LogicalVolume* lTile = new G4LogicalVolume(tileBox, tileMaterial, "Tile", 0, ilayer == 0 ? m_sensitiveAero : 0);
581  setColor(*lTile, "rgb(0.0, 1.0, 1.0,1.0)");
582  G4Transform3D trans = G4Translate3D(posX, posY, posZ + sizeZ / 2. - zFrame / 2.);
583  new G4PVPlacement(trans, lTile, "ARICH.tile", lenvBox, false, ilayer + 1);
584  posZ += sizeZ;
585  }
586  if (m_arichbtgp->getAverageAgel() && m_agelthickness.size()) {
587  B2INFO("Average aerogel will be used in the reconstruction ");
588  m_arichgp->setAeroRefIndex(0, meanrefind / m_agelthickness.size());
589  m_arichgp->setAerogelZPosition(0, (posZ0 - zFrame)* Unit::mm / CLHEP::mm);
590  m_arichgp->setAerogelThickness(0, posZ * Unit::mm / CLHEP::mm);
591  if (meantrlen > 0 && posZ > 0) meantrlen = 1 / meantrlen / posZ;
592  m_arichgp->setAeroTransLength(0, meantrlen);
593  }
594 
595 
596  // place mirrors
597  GearDir mirrorsParam(mirrorcontent, "Mirrors");
598  double height = mirrorsParam.getLength("height") * CLHEP::mm / Unit::mm;
599  double width = mirrorsParam.getLength("width") * CLHEP::mm / Unit::mm;
600  double thickness = mirrorsParam.getLength("thickness") * CLHEP::mm / Unit::mm;
601  string mirrMat = mirrorsParam.getString("material");
602  G4Material* mirrMaterial = Materials::get(mirrMat);
603  G4Box* mirrBox = new G4Box("mirrBox", thickness / 2., height / 2., width / 2.);
604  G4LogicalVolume* lmirror = new G4LogicalVolume(mirrBox, mirrMaterial, "mirror");
605 
606  Materials& materials = Materials::getInstance();
607  GearDir surface(mirrorsParam, "Surface");
608  G4OpticalSurface* optSurf = materials.createOpticalSurface(surface);
609  new G4LogicalSkinSurface("mirrorsSurface", lmirror, optSurf);
610  int iMirror = 0;
611  BOOST_FOREACH(const GearDir & mirror, mirrorsParam.getNodes("Mirror")) {
612  double xpos = mirror.getLength("xPos") * CLHEP::mm / Unit::mm;
613  double ypos = mirror.getLength("yPos") * CLHEP::mm / Unit::mm;
614  double zpos = mirror.getLength("zPos") * CLHEP::mm / Unit::mm;
615  double angle = mirror.getAngle("angle") / Unit::rad;
616  G4ThreeVector origin(xpos, ypos, zpos + width / 2. - zFrame / 2.);
617  G4RotationMatrix Ra;
618  Ra.rotateZ(angle);
619  G4Transform3D trans = G4Transform3D(Ra, origin);
620  new G4PVPlacement(G4Transform3D(Ra, origin), lmirror, "ARICH.mirror", lenvBox, false, iMirror);
621  iMirror++;
622  }
623  m_arichgp->Print();
624  m_arichbtgp->Print();
625  }
626 
627 
628  }
630 }
Belle2::ARICHBtestGeometryPar
The Class for ARICH Beamtest Geometry Parameters.
Definition: ARICHBtestGeometryPar.h:99
Belle2::ARICHGeometryPar
The Class for ARICH Geometry Parameters.
Definition: ARICHGeometryPar.h:35
Belle2::ARICHTracking
Beamtest ARICH Geometry Tracking Class.
Definition: ARICHBtestGeometryPar.h:26
Belle2::ARICHTracking::slp
float slp[2]
Calibration constants of the MWPC (\delta x= slope \delta t + offset) - slopes for x an y direction.
Definition: ARICHBtestGeometryPar.h:54
Belle2::ARICHTracking::cutll
int cutll[2]
Cuts on the tdc sums - lower levels.
Definition: ARICHBtestGeometryPar.h:58
Belle2::ARICHTracking::pos
float pos[3]
MWPC chamber position.
Definition: ARICHBtestGeometryPar.h:62
Belle2::ARICHTracking::tdc
int tdc[4]
TDC of the 4 cathode signals.
Definition: ARICHBtestGeometryPar.h:48
Belle2::ARICHTracking::atdc
int atdc
TDC of the anode signal.
Definition: ARICHBtestGeometryPar.h:52
Belle2::ARICHTracking::offset
float offset[2]
Calibration constants of the MWPC - offsets for x an y direction.
Definition: ARICHBtestGeometryPar.h:56
Belle2::ARICHTracking::cutul
int cutul[2]
Cuts on the tdc sums - upper levels.
Definition: ARICHBtestGeometryPar.h:60
Belle2::GearDir
GearDir is the basic class used for accessing the parameter store.
Definition: GearDir.h:31
Belle2::GearDir::getString
virtual std::string getString(const std::string &path="") const noexcept(false) override
Get the parameter path as a string.
Definition: GearDir.h:69
Belle2::Module
Base class for Modules.
Definition: Module.h:72
Belle2::StoreObjPtr
Type-safe access to single objects in the data store.
Definition: StoreObjPtr.h:96
Belle2::Unit::mm
static const double mm
[millimeters]
Definition: Unit.h:70
Belle2::Unit::rad
static const double rad
Standard of [angle].
Definition: Unit.h:50
Belle2::Unit::eV
static const double eV
[electronvolt]
Definition: Unit.h:112
Belle2::Unit::MeV
static const double MeV
[megaelectronvolt]
Definition: Unit.h:114
Belle2::arich::GeoARICHBtestCreator::m_comment
std::string m_comment
comment in the runlog
Definition: GeoARICHBtestCreator.h:100
Belle2::arich::GeoARICHBtestCreator::m_aerosupport
int m_aerosupport
Type of aerogel support - not used at the moment.
Definition: GeoARICHBtestCreator.h:105
Belle2::arich::GeoARICHBtestCreator::m_aerogelID
std::string m_aerogelID
ID of the aerogel configuration setup.
Definition: GeoARICHBtestCreator.h:86
Belle2::arich::GeoARICHBtestCreator::m_rotation1
double m_rotation1
rotation angle of the frame
Definition: GeoARICHBtestCreator.h:111
Belle2::arich::GeoARICHBtestCreator::m_runtype
std::string m_runtype
Type of the beamtest run.
Definition: GeoARICHBtestCreator.h:82
Belle2::arich::GeoARICHBtestCreator::m_agelrefind
std::vector< double > m_agelrefind
vector of aerogel refractive indices
Definition: GeoARICHBtestCreator.h:118
Belle2::arich::GeoARICHBtestCreator::createBtestGeometry
void createBtestGeometry(const GearDir &content, G4LogicalVolume &topVolume)
Creation of the beamtest geometry.
Definition: GeoARICHBtestCreator.cc:356
Belle2::arich::GeoARICHBtestCreator::m_mytype
std::string m_mytype
type of the run
Definition: GeoARICHBtestCreator.h:96
Belle2::arich::GeoARICHBtestCreator::createAerogel
G4Material * createAerogel(const char *aeroname, double rind, double trl)
create aerogel material
Definition: GeoARICHBtestCreator.cc:302
Belle2::arich::GeoARICHBtestCreator::m_ry
double m_ry
y shift of the prototype ARICH frame
Definition: GeoARICHBtestCreator.h:94
Belle2::arich::GeoARICHBtestCreator::m_ageltrlen
std::vector< double > m_ageltrlen
vector of aerogel transmission lengths
Definition: GeoARICHBtestCreator.h:122
Belle2::arich::GeoARICHBtestCreator::~GeoARICHBtestCreator
virtual ~GeoARICHBtestCreator()
The destructor of the GeoPXDCreator class.
Definition: GeoARICHBtestCreator.cc:83
Belle2::arich::GeoARICHBtestCreator::buildModule
G4LogicalVolume * buildModule(GearDir Module)
Build the module.
Definition: GeoARICHBtestCreator.cc:221
Belle2::arich::GeoARICHBtestCreator::getAvgRINDEX
double getAvgRINDEX(G4Material *material)
Get the average refractive index if the material.
Definition: GeoARICHBtestCreator.cc:211
Belle2::arich::GeoARICHBtestCreator::m_mirrorID
std::string m_mirrorID
ID of the mirror configuration setup.
Definition: GeoARICHBtestCreator.h:88
Belle2::arich::GeoARICHBtestCreator::create
virtual void create(const GearDir &content, G4LogicalVolume &topVolume, geometry::GeometryTypes type)
Creates the ROOT Objects for the ARICH Beamtest 2011 geometry.
Definition: GeoARICHBtestCreator.cc:89
Belle2::arich::GeoARICHBtestCreator::m_aerogeldx
double m_aerogeldx
shift of the aerogel center
Definition: GeoARICHBtestCreator.h:107
Belle2::arich::GeoARICHBtestCreator::m_datum
std::string m_datum
datum of the runlog
Definition: GeoARICHBtestCreator.h:102
Belle2::arich::GeoARICHBtestCreator::m_sensitiveAero
SensitiveAero * m_sensitiveAero
pointer to the sesnitive aerogel
Definition: GeoARICHBtestCreator.h:69
Belle2::arich::GeoARICHBtestCreator::m_author
std::string m_author
Beamtest runlog record author.
Definition: GeoARICHBtestCreator.h:78
Belle2::arich::GeoARICHBtestCreator::m_daqqa
std::string m_daqqa
classification of the run
Definition: GeoARICHBtestCreator.h:98
Belle2::arich::GeoARICHBtestCreator::m_hapdID
std::string m_hapdID
ID of the HAPD configuration setup.
Definition: GeoARICHBtestCreator.h:84
Belle2::arich::GeoARICHBtestCreator::m_agelthickness
std::vector< double > m_agelthickness
vector of aerogel thicknesses
Definition: GeoARICHBtestCreator.h:120
Belle2::arich::GeoARICHBtestCreator::m_rotation
double m_rotation
rotation angle of the setup
Definition: GeoARICHBtestCreator.h:90
Belle2::arich::GeoARICHBtestCreator::m_sensitive
SensitiveDetector * m_sensitive
pointer to the sensitive detector
Definition: GeoARICHBtestCreator.h:67
Belle2::arich::GeoARICHBtestCreator::m_rx
double m_rx
x shift of the prototype ARICH frame
Definition: GeoARICHBtestCreator.h:92
Belle2::arich::GeoARICHBtestCreator::m_neve
int m_neve
Number of event in the beamtest run.
Definition: GeoARICHBtestCreator.h:80
Belle2::arich::GeoARICHBtestCreator::m_configuration
int m_configuration
configuration number of the HAPD
Definition: GeoARICHBtestCreator.h:113
Belle2::arich::SensitiveAero
This is optional (temporary) class that provides information on track parameters on aerogel plane,...
Definition: SensitiveAero.h:22
Belle2::arich::SensitiveDetector
The Class for ARICH Sensitive Detector.
Definition: SensitiveDetector.h:26
Belle2::gearbox::Interface::getAngle
double getAngle(const std::string &path="") const noexcept(false)
Get the parameter path as a double converted to the standard angle unit.
Definition: Interface.h:299
Belle2::gearbox::Interface::getArray
std::vector< double > getArray(const std::string &path) const noexcept(false)
Get the parameter path as a list of double values converted to the standard unit.
Definition: Interface.cc:123
Belle2::gearbox::Interface::getDouble
double getDouble(const std::string &path="") const noexcept(false)
Get the parameter path as a double.
Definition: Interface.cc:41
Belle2::gearbox::Interface::getPath
std::string getPath() const
Return path of the current interface.
Definition: Interface.h:70
Belle2::gearbox::Interface::getLength
double getLength(const std::string &path="") const noexcept(false)
Get the parameter path as a double converted to the standard length unit.
Definition: Interface.h:259
Belle2::gearbox::Interface::getNodes
std::vector< GearDir > getNodes(const std::string &path="") const
Get vector of GearDirs which point to all the nodes the given path evaluates to.
Definition: Interface.cc:21
Belle2::gearbox::Interface::getInt
int getInt(const std::string &path="") const noexcept(false)
Get the parameter path as a int.
Definition: Interface.cc:60
Belle2::geometry::Materials
Thin wrapper around the Geant4 Material system.
Definition: Materials.h:48
Belle2::geometry::Materials::get
static G4Material * get(const std::string &name)
Find given material.
Definition: Materials.h:63
Belle2::geometry::Materials::getInstance
static Materials & getInstance()
Get a reference to the singleton instance.
Definition: Materials.cc:85
Belle2::geometry::Materials::createOpticalSurface
G4OpticalSurface * createOpticalSurface(const gearbox::Interface &parameters)
Create an optical surface from parameters, will abort on error.
Definition: Materials.cc:291
Belle2::geometry::Materials::getElement
G4Element * getElement(const std::string &name)
Find given chemical element.
Definition: Materials.cc:150
LogVar
Class to store variables with their name which were sent to the logging service.
Definition: LogVariableStream.h:24
Belle2::ARICHGeometryPar::getDetectorXPadNumber
int getDetectorXPadNumber()
get number of pads of detector module (in one direction)
Definition: ARICHGeometryPar.h:355
Belle2::ARICHBtestGeometryPar::setRotationCenter
void setRotationCenter(const TVector3 &)
Set the rotation center of the Aerogel RICH frame.
Definition: ARICHBtestGeometryPar.cc:117
Belle2::ARICHGeometryPar::setAeroTransLength
void setAeroTransLength(int ilayer, double trlen)
set transmission length of "ilayer" aerogel layer
Definition: ARICHGeometryPar.cc:511
Belle2::ARICHBtestGeometryPar::AddHapdElectronicMapPair
int AddHapdElectronicMapPair(int, int)
Set the mapping of the electronic channel to the HAPD module nr and the channel number.
Definition: ARICHBtestGeometryPar.cc:97
Belle2::ARICHGeometryPar::setAerogelThickness
void setAerogelThickness(int ilayer, double thick)
set thickness of "ilayer" aerogel layer
Definition: ARICHGeometryPar.cc:521
Belle2::ARICHGeometryPar::getSensitiveSurfaceSize
double getSensitiveSurfaceSize() const
get size of detector sensitive surface (size of two chips + gap between)
Definition: ARICHGeometryPar.h:284
Belle2::ARICHBtestGeometryPar::setTrackingShift
void setTrackingShift(const TVector3 &)
Set the tracking shift.
Definition: ARICHBtestGeometryPar.cc:103
Belle2::ARICHBtestGeometryPar::getAverageAgel
bool getAverageAgel()
Get the flag for the reconstruction by using the average aerogel refractive index.
Definition: ARICHBtestGeometryPar.cc:148
Belle2::ARICHBtestGeometryPar::setFrameRotation
void setFrameRotation(double)
Set the rotation angle of the Aerogel RICH frame.
Definition: ARICHBtestGeometryPar.cc:130
Belle2::ARICHBtestGeometryPar::setMwpc
void setMwpc(ARICHTracking *m_mwpc)
Set the pointer of the tracking MWPCs.
Definition: ARICHBtestGeometryPar.cc:70
Belle2::ARICHGeometryPar::Initialize
void Initialize(const GearDir &content)
calculates detector parameters needed for geometry build and reconstruction.
Definition: ARICHGeometryPar.cc:70
Belle2::ARICHGeometryPar::setAerogelZPosition
void setAerogelZPosition(int ilayer, double zPos)
set z position of "ilayer" aerogel layer
Definition: ARICHGeometryPar.cc:527
Belle2::ARICHBtestGeometryPar::Instance
static ARICHBtestGeometryPar * Instance()
Static method to get a reference to the ARICHBtestGeometryPar instance.
Definition: ARICHBtestGeometryPar.cc:27
Belle2::ARICHGeometryPar::Print
void Print(void) const
Print some debug information.
Definition: ARICHGeometryPar.cc:302
Belle2::ARICHGeometryPar::Instance
static ARICHGeometryPar * Instance()
Static method to get a reference to the ARICHGeometryPar instance.
Definition: ARICHGeometryPar.cc:32
Belle2::ARICHGeometryPar::setAeroRefIndex
void setAeroRefIndex(int ilayer, double n)
set refractive index of "ilayer" aerogel layer
Definition: ARICHGeometryPar.cc:516
Belle2::ARICHGeometryPar::getModAngle
double getModAngle(int copyno)
get the angle of copyno-th HAPD rotation
Definition: ARICHGeometryPar.cc:410
Belle2::ARICHGeometryPar::getOriginG4
G4ThreeVector getOriginG4(int copyNo)
get the position of copyNo-th HAPD module origin (returns G4ThreeVector)
Definition: ARICHGeometryPar.cc:404
Belle2::ARICHBtestGeometryPar::setOffset
void setOffset(const TVector3 &)
Set of the setup global offset.
Definition: ARICHBtestGeometryPar.cc:164
Belle2::ARICHBtestGeometryPar::setAverageAgel
void setAverageAgel(bool)
Set the flag for the reconstruction by using the average aerogel refractive index.
Definition: ARICHBtestGeometryPar.cc:143
Belle2::ARICHGeometryPar::getNMCopies
int getNMCopies() const
get the total number of HAPD modules
Definition: ARICHGeometryPar.h:289
Belle2::ARICHGeometryPar::setActive
void setActive(int module, int channel, bool val)
set the channel on/off
Definition: ARICHGeometryPar.cc:547
Belle2::ARICHBtestGeometryPar::AddHapdChannelPositionPair
int AddHapdChannelPositionPair(double, double)
Set the position of the HAPD channel.
Definition: ARICHBtestGeometryPar.cc:91
Belle2::ARICHGeometryPar::setWindowRefIndex
void setWindowRefIndex(double refInd)
set detector module window refractive index
Definition: ARICHGeometryPar.cc:532
Belle2::geometry::setVisibility
void setVisibility(G4LogicalVolume &volume, bool visible)
Helper function to quickly set the visibility of a given volume.
Definition: utilities.cc:108
Belle2::geometry::setColor
void setColor(G4LogicalVolume &volume, const std::string &color)
Set the color of a logical volume.
Definition: utilities.cc:100
Belle2::geometry::GeometryTypes
GeometryTypes
Flag indiciating the type of geometry to be used.
Definition: GeometryManager.h:36
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:17