From d302838e0a8573fab5b0bdd9a535a61f0d003708 Mon Sep 17 00:00:00 2001
From: David Emschermann <d.emschermann@gsi.de>
Date: Wed, 9 Mar 2022 09:52:33 +0100
Subject: [PATCH] add RICH v21b, which is shifted 20 cm along the rails in -x
directions, refs #2459
---
.../geometry/rich/create_rich_v21b_mcbm.C | 370 ++++++++++++++++++
macro/mcbm/mcbm_transport.C | 5 +-
2 files changed, 373 insertions(+), 2 deletions(-)
create mode 100644 macro/mcbm/geometry/rich/create_rich_v21b_mcbm.C
diff --git a/macro/mcbm/geometry/rich/create_rich_v21b_mcbm.C b/macro/mcbm/geometry/rich/create_rich_v21b_mcbm.C
new file mode 100644
index 0000000000..d0dff0d0e2
--- /dev/null
+++ b/macro/mcbm/geometry/rich/create_rich_v21b_mcbm.C
@@ -0,0 +1,370 @@
+/* Copyright (C) 2020-2022 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt
+ SPDX-License-Identifier: GPL-3.0-only
+ Authors: Florian Uhlig, David Emschermann [committer] */
+
+#include <iostream>
+//#include "FairGeoMedium.h"
+//#include "FairGeoBuilder.h"
+//#include "FairGeoMedia.h"
+#include "TGeoManager.h"
+#include "TGeoVolume.h"
+
+using namespace std;
+
+// Changelog
+
+// 2022-03-08 - v21b - DE - Based on v21a, relocate RICH along the rail onto the 25 degree line
+// 2021-09-28 - v21a - SR - Based on v20d, moved RICH according to the TOF v21d
+// 2020-05-25 - v20d - AW - Move mRICH by +7.2cm in Z accordingly to the mTOF movement from rev 16441
+// 2020-05-20 - v20d - AW - changed positioning to +2.5cm to 25°line (nearer to beam axis) (from Analysis)
+// 2020-05-06 - v20b - AW - changed positioning to center on 25° line (without Window)
+// 2020-04-02 - v19a - FU - same as v19c but with exported geometry in output file
+// 2019-11-28 - v19c - DE - move mRICH +12+7 cm in x direction (same as mTOF v19b) for the Nov 2019 run
+// 2017-11-17 - v18d - DE - add aerogel as radiator to the mRICH module
+// 2017-11-12 - v18c - DE - push mRICH downstream to z=355 cm for long setup
+// 2017-07-19 - v18b - DE - add one level to the geometry hierarchy
+// 2017-07-15 - v18b - DE - arrange 4x mRICH detectors in the setup and tilt towards target
+
+
+void create_rich_v21b_mcbm()
+{
+ gSystem->Load("libGeom");
+
+ TString geoFileName = "rich_v21b_mcbm.geo.root";
+ TString topNodeName = "rich_v21b_mcbm";
+
+ Double_t XOffset = 0;
+ Double_t ZOffset = 0;
+
+ // displace RICH module along the rails on the TOF frame
+ Double_t displacement_along_rails = -20.0; // cm for Carbon-12 run
+ Double_t support_frame_rotation_angle = 12.5; // degrees
+
+ XOffset = +displacement_along_rails *cos( support_frame_rotation_angle * acos(-1)/180 );
+ ZOffset = -displacement_along_rails *sin( support_frame_rotation_angle * acos(-1)/180 );
+ // end of displacement along the rails
+
+ FairGeoLoader* geoLoad = new FairGeoLoader("TGeo", "FairGeoLoader");
+ FairGeoInterface* geoFace = geoLoad->getGeoInterface();
+ TString geoPath = gSystem->Getenv("VMCWORKDIR");
+ TString medFile = geoPath + "/geometry/media.geo";
+
+ geoFace->setMediaFile(medFile);
+ geoFace->readMedia();
+ TGeoManager* gGeoMan = gGeoManager;
+
+ FairGeoMedia* geoMedia = geoFace->getMedia();
+ FairGeoBuilder* geoBuild = geoLoad->getGeoBuilder();
+
+ //Media
+ FairGeoMedium* mAluminium = geoMedia->getMedium("aluminium");
+ if (mAluminium == NULL) Fatal("Main", "FairMedium aluminium not found");
+ geoBuild->createMedium(mAluminium);
+ TGeoMedium* medAl = gGeoMan->GetMedium("aluminium");
+ if (medAl == NULL) Fatal("Main", "Medium vacuum not found");
+
+ FairGeoMedium* mKapton = geoMedia->getMedium("kapton");
+ if (mKapton == NULL) Fatal("Main", "FairMedium kapton not found");
+ geoBuild->createMedium(mKapton);
+ TGeoMedium* medKapton = gGeoMan->GetMedium("kapton");
+ if (medKapton == NULL) Fatal("Main", "Medium kapton not found");
+
+ FairGeoMedium* mAir = geoMedia->getMedium("air");
+ if (mAir == NULL) Fatal("Main", "FairMedium air not found");
+ geoBuild->createMedium(mAir);
+ TGeoMedium* medAir = gGeoMan->GetMedium("air");
+ if (medAir == NULL) Fatal("Main", "Medium air not found");
+
+ FairGeoMedium* mPmtGlas = geoMedia->getMedium("PMTglass");
+ if (mPmtGlas == NULL) Fatal("Main", "FairMedium vacuum not found");
+ geoBuild->createMedium(mPmtGlas);
+ TGeoMedium* medPmtGlas = gGeoMan->GetMedium("PMTglass");
+ if (medPmtGlas == NULL) Fatal("Main", "Medium vacuum not found");
+
+ FairGeoMedium* mNitrogen = geoMedia->getMedium("RICHgas_N2_dis");
+ if (mNitrogen == NULL) Fatal("Main", "FairMedium RICHgas_N2_dis not found");
+ geoBuild->createMedium(mNitrogen);
+ TGeoMedium* medNitrogen = gGeoMan->GetMedium("RICHgas_N2_dis");
+ if (medNitrogen == NULL) Fatal("Main", "Medium RICHgas_N2_dis not found");
+
+ FairGeoMedium* mPMT = geoMedia->getMedium("CsI");
+ if (mPMT == NULL) Fatal("Main", "FairMedium CsI not found");
+ geoBuild->createMedium(mPMT);
+ TGeoMedium* medCsI = gGeoMan->GetMedium("CsI");
+ if (medCsI == NULL) Fatal("Main", "Medium CsI not found");
+
+ FairGeoMedium* mElectronic = geoMedia->getMedium("air");
+ if (mElectronic == NULL) Fatal("Main", "FairMedium air not found");
+ geoBuild->createMedium(mElectronic);
+ TGeoMedium* medElectronic = gGeoMan->GetMedium("air");
+ if (medElectronic == NULL) Fatal("Main", "Medium air not found");
+
+ FairGeoMedium* mAerogel = geoMedia->getMedium("aerogel");
+ if (mAerogel == NULL) Fatal("Main", "FairMedium aerogel not found");
+ geoBuild->createMedium(mAerogel);
+ TGeoMedium* medAerogel = gGeoMan->GetMedium("aerogel");
+ if (medAerogel == NULL) Fatal("Main", "Medium aerogel not found");
+
+
+ //Dimensions of the RICH prototype [cm]
+
+ const Double_t ItemToKapton = 1.85;
+
+ // Box
+ const Double_t boxThickness = 0.3;
+ const Double_t boxBackThickness = 1.;
+ const Double_t boxFrontThickness = 0.3;
+ const Double_t boxXLen = 38.2;
+ const Double_t boxYLen = 59.7;
+ const Double_t boxZLen = 25. + boxBackThickness - ItemToKapton;
+ const Bool_t isRotate90DegZ = false;
+
+ // PMT
+ // PMT specs https://www.hamamatsu.com/resources/pdf/etd/H12700_TPMH1348E.pdf
+ const Double_t pmtNofPixels = 8;
+ const Double_t pmtSize = 5.2;
+ const Double_t pmtSizeHalf = pmtSize / 2.;
+ const Double_t pmtSensitiveSize = 4.85;
+ const Double_t pmtPixelSize = 0.6; //pmtSensitiveSize / pmtNofPixels;
+ const Double_t pmtPixelSizeHalf = pmtPixelSize / 2.;
+ const Double_t pmtOuterPixelSize = 0.625; //pmtSensitiveSize / pmtNofPixels;
+ const Double_t pmtOuterPixelSizeHalf = pmtOuterPixelSize / 2.;
+ const Double_t pmtGap = 0.1;
+ const Double_t pmtGapHalf = pmtGap / 2.;
+ const Double_t pmtMatrixGap = 0.1;
+ const Double_t pmtPixelThickness = 0.1;
+ const Double_t pmtWindowThickness = 0.00; //0.15
+ Double_t pmtWindowPixelDistance = 0.0;
+ if (pmtWindowThickness > 0.0) pmtWindowPixelDistance = 0.0;
+ const Double_t pmtThickness = pmtPixelThickness + pmtWindowThickness + pmtWindowPixelDistance;
+ const Double_t pmtZPadding = 3.871;
+ const Double_t pmtXShift = 2.5;
+
+ // Electronics
+ const Double_t elecLength = 10.;
+ const Double_t elecWidth = 3. * pmtSize;
+ const Double_t elecHeight = 2. * pmtSize;
+
+ // Aerogel Box
+ const Double_t aerogelXLen = 20.;
+ const Double_t aerogelYLen = 20.;
+ const Double_t aerogelZLen = 3.;
+ const Double_t aerogelPmtDistance = 10.;
+ const Double_t aerogelGap = 0.1;
+
+ // Imaginary sensitive plane
+ Bool_t isIncludeSensPlane = false;
+ const Double_t sensPlaneSize = 200.;
+ const Double_t sensPlaneBoxDistance = 1.;
+
+ //Global positioning
+ const Double_t boxFrontToTarget = 334.2 + ItemToKapton + boxZLen / 2;
+
+ // Create the top volume
+ TGeoVolume* top = new TGeoVolumeAssembly("top");
+ gGeoMan->SetTopVolume(top);
+
+ // Rich assembly
+ TGeoVolume* rich = new TGeoVolumeAssembly(topNodeName);
+ top->AddNode(rich, 1);
+
+ // Al box
+ TGeoRotation* rotBox = new TGeoRotation("rotBox", 0., 0., 0.);
+ TGeoCombiTrans* trBox = NULL;
+ if (isRotate90DegZ) {
+ rotBox->RotateZ(90.);
+ rotBox->RotateY(0.955); //Tilting of Rich
+ trBox = new TGeoCombiTrans(0, -pmtXShift - (pmtSize + pmtGap + 1), boxFrontToTarget, rotBox);
+ }
+ else {
+ rotBox->RotateX(0.955); //Tilting of Rich
+ // Nov 2019 run
+ rotBox->RotateY(support_frame_rotation_angle);
+ //trBox = new TGeoCombiTrans(-pmtXShift-(pmtSize + pmtGap + 1)+12.+7., 0, boxFrontToTarget, rotBox);
+ // Mar2020
+ trBox = new TGeoCombiTrans(-pmtXShift - (pmtSize + pmtGap + 1) + 6.3 + XOffset + 16. + 5. + 2.5, 0,
+ boxFrontToTarget + 27. + ZOffset, rotBox);
+ }
+ TGeoVolume* boxVol = gGeoMan->MakeBox("box", medAl, boxXLen / 2., boxYLen / 2., boxZLen / 2.);
+ // rich->AddNode(boxVol, 1, trBox);
+ rich->AddNode(boxVol, 1);
+
+ // Gas
+ TGeoVolume* gasVol = gGeoMan->MakeBox("Gas", medNitrogen, boxXLen / 2. - boxThickness, boxYLen / 2. - boxThickness,
+ ((boxZLen - boxBackThickness - boxFrontThickness) / 2.));
+ TGeoTranslation* trGas = new TGeoTranslation(0., 0., (-(boxBackThickness - boxFrontThickness) / 2.));
+ boxVol->AddNode(gasVol, 1, trGas);
+
+ // Front plane from kapton
+ TGeoVolume* kaptonVol = gGeoMan->MakeBox("kapton", medKapton, boxXLen / 2. - boxThickness,
+ boxYLen / 2. - boxThickness, boxFrontThickness / 2.);
+ TGeoTranslation* trKapton = new TGeoTranslation(0., 0., -(boxZLen - boxFrontThickness) / 2.);
+ boxVol->AddNode(kaptonVol, 1, trKapton);
+
+
+ // Aerogel
+ TGeoVolume* aerogelVol =
+ gGeoMan->MakeBox("aerogel", medAerogel, aerogelXLen / 2., aerogelYLen / 2., aerogelZLen / 2.);
+ Double_t aerogelZ = boxZLen / 2. - boxBackThickness - pmtZPadding - aerogelPmtDistance - aerogelZLen / 2.;
+ TGeoTranslation* trAerogel1 = new TGeoTranslation(pmtXShift, (aerogelYLen + aerogelGap) / 2., aerogelZ);
+ TGeoTranslation* trAerogel2 = new TGeoTranslation(pmtXShift, -(aerogelYLen + aerogelGap) / 2., aerogelZ);
+ gasVol->AddNode(aerogelVol, 1, trAerogel1);
+ gasVol->AddNode(aerogelVol, 2, trAerogel2);
+
+ // PMT
+ TGeoTranslation* trPmt1 = new TGeoTranslation(pmtSizeHalf + pmtGapHalf, -pmtSize - pmtGap, 0.); // Bottom Right
+ TGeoTranslation* trPmt2 = new TGeoTranslation(pmtSizeHalf + pmtGapHalf, 0., 0.); // Middle Right
+ TGeoTranslation* trPmt3 = new TGeoTranslation(pmtSizeHalf + pmtGapHalf, pmtSize + pmtGap, 0.); // Top Right
+ TGeoTranslation* trPmt4 = new TGeoTranslation(-pmtSizeHalf - pmtGapHalf, -pmtSize - pmtGap, 0.); // Bottom Left
+ TGeoTranslation* trPmt5 = new TGeoTranslation(-pmtSizeHalf - pmtGapHalf, 0., 0.); // Middle Left
+ TGeoTranslation* trPmt6 = new TGeoTranslation(-pmtSizeHalf - pmtGapHalf, pmtSize + pmtGap, 0.); // Top Left
+
+ TGeoVolume* pmtContVol = new TGeoVolumeAssembly("pmt_cont_vol");
+ TGeoVolume* pmtVol = new TGeoVolumeAssembly("pmt_vol_0");
+ TGeoVolume* pmtVol_180 = new TGeoVolumeAssembly("pmt_vol_1");
+ TGeoVolume* pmtPixelVol =
+ gGeoMan->MakeBox("pmt_pixel", medCsI, pmtPixelSize / 2., pmtPixelSize / 2., pmtPixelThickness / 2.);
+ TGeoVolume* pmtPixelVol_outer_outer =
+ gGeoMan->MakeBox("pmt_pixel", medCsI, pmtOuterPixelSize / 2., pmtOuterPixelSize / 2., pmtPixelThickness / 2.);
+ TGeoVolume* pmtPixelVol_inner_outer =
+ gGeoMan->MakeBox("pmt_pixel", medCsI, pmtPixelSize / 2., pmtOuterPixelSize / 2., pmtPixelThickness / 2.);
+ TGeoVolume* pmtPixelVol_outer_inner =
+ gGeoMan->MakeBox("pmt_pixel", medCsI, pmtOuterPixelSize / 2., pmtPixelSize / 2., pmtPixelThickness / 2.);
+ TGeoVolume* pmtWindow = gGeoMan->MakeBox("pmt_Window", medPmtGlas, pmtSizeHalf, pmtSizeHalf, pmtWindowThickness / 2.);
+
+ pmtContVol->AddNode(pmtVol_180, 1, trPmt1); // To Power Module // Bottom Right
+ pmtContVol->AddNode(pmtVol, 2, trPmt2); // Middle // Middle Right
+ pmtContVol->AddNode(pmtVol, 3, trPmt3); // To Combiner // Top Right
+ pmtContVol->AddNode(pmtVol_180, 4, trPmt4); // To Power Module // Bottom Left
+ pmtContVol->AddNode(pmtVol, 5, trPmt5); // Middle // Middle Left
+ pmtContVol->AddNode(pmtVol, 6, trPmt6); // To Combiner // Top Left
+
+ Int_t halfPmtNofPixels = pmtNofPixels / 2;
+ for (Int_t i = 0; i < pmtNofPixels; i++) {
+ for (Int_t j = 0; j < pmtNofPixels; j++) {
+ //check outer Pixels
+
+ Double_t x = 0.;
+ Double_t y = 0.;
+ bool x_outer = false;
+ bool y_outer = false;
+
+ if (i == 0) {
+ x = -(halfPmtNofPixels - 1) * pmtPixelSize - pmtOuterPixelSizeHalf;
+ x_outer = true;
+ }
+ else if (i == (pmtNofPixels - 1)) {
+ x = (halfPmtNofPixels - 1) * pmtPixelSize + pmtOuterPixelSizeHalf;
+ x_outer = true;
+ }
+ else if (i < halfPmtNofPixels) {
+ x = -((halfPmtNofPixels - 1) - i) * pmtPixelSize - pmtPixelSizeHalf;
+ x_outer = false;
+ }
+ else {
+ x = (i - halfPmtNofPixels) * pmtPixelSize + pmtPixelSizeHalf;
+ x_outer = false;
+ }
+
+
+ if (j == 0) {
+ y = -(halfPmtNofPixels - 1) * pmtPixelSize - pmtOuterPixelSizeHalf;
+ y_outer = true;
+ }
+ else if (j == (pmtNofPixels - 1)) {
+ y = (halfPmtNofPixels - 1) * pmtPixelSize + pmtOuterPixelSizeHalf;
+ y_outer = true;
+ }
+ else if (j < halfPmtNofPixels) {
+ y = -((halfPmtNofPixels - 1) - j) * pmtPixelSize - pmtPixelSizeHalf;
+ y_outer = false;
+ }
+ else {
+ y = (j - halfPmtNofPixels) * pmtPixelSize + pmtPixelSizeHalf;
+ y_outer = false;
+ }
+ //Double_t x = (i - (halfPmtNofPixels - 1)) * pmtPixelSize - pmtPixelSizeHalf;
+
+ //Double_t y = (j - (halfPmtNofPixels - 1)) * pmtPixelSize - pmtPixelSizeHalf;
+ TGeoTranslation* trij = new TGeoTranslation(x, y, -pmtPixelThickness / 2. + pmtThickness / 2.);
+
+ //normal rotated PMT (HV connector at lower side)
+ uint uid = 18 - ((j % 4) / 2) * 16 - (j % 2) + 2 * i + (1 - (j / 4)) * 100;
+
+ //180?? rotated PMT
+ uint uid_180 = 15 + ((j % 4) / 2) * 16 + (j % 2) - 2 * i + (1 - (j / 4)) * 100;
+ if (x_outer == true && y_outer == true) {
+ pmtVol->AddNode(pmtPixelVol_outer_outer, uid, trij);
+ pmtVol_180->AddNode(pmtPixelVol_outer_outer, uid_180, trij);
+ }
+ else if (x_outer == true && y_outer == false) {
+ pmtVol->AddNode(pmtPixelVol_outer_inner, uid, trij);
+ pmtVol_180->AddNode(pmtPixelVol_outer_inner, uid_180, trij);
+ }
+ else if (x_outer == false && y_outer == true) {
+ pmtVol->AddNode(pmtPixelVol_inner_outer, uid, trij);
+ pmtVol_180->AddNode(pmtPixelVol_inner_outer, uid_180, trij);
+ }
+ else {
+ pmtVol->AddNode(pmtPixelVol, uid, trij);
+ pmtVol_180->AddNode(pmtPixelVol, uid_180, trij);
+ }
+ }
+ }
+
+ TGeoTranslation* trWndw = new TGeoTranslation(0., 0., +pmtWindowThickness / 2. - pmtThickness / 2.);
+ if (pmtWindowThickness > 0.0) {
+ pmtVol->AddNode(pmtWindow, 1, trWndw);
+ pmtVol_180->AddNode(pmtWindow, 1, trWndw);
+ }
+
+ //
+ Double_t pmtContXLen = 2 * (pmtSize + pmtGap);
+ Double_t pmtContYLen = 3 * (pmtSize + pmtGap);
+ TGeoTranslation* trPmtCont1 = new TGeoTranslation(-pmtContXLen / 2., pmtContYLen, 0.); //Top Left
+ TGeoTranslation* trPmtCont2 = new TGeoTranslation(-pmtContXLen / 2., 0., 0.); //Middle Left
+ TGeoTranslation* trPmtCont3 = new TGeoTranslation(-pmtContXLen / 2., -pmtContYLen, 0.); //Bottom Left
+ TGeoTranslation* trPmtCont4 = new TGeoTranslation(pmtContXLen / 2., pmtContYLen, 0.); //
+ TGeoTranslation* trPmtCont5 = new TGeoTranslation(pmtContXLen / 2., 0., 0.); //
+ TGeoTranslation* trPmtCont6 = new TGeoTranslation(pmtContXLen / 2., -pmtContYLen, 0.); //
+
+ TGeoVolume* pmtPlaneVol = new TGeoVolumeAssembly("pmt_plane");
+ pmtPlaneVol->AddNode(pmtContVol, 00,
+ trPmtCont1); // Id gives the Address 0x7**0 of the top left DiRICH in a BP (X|Y)
+ pmtPlaneVol->AddNode(pmtContVol, 03, trPmtCont2);
+ pmtPlaneVol->AddNode(pmtContVol, 06, trPmtCont3);
+ pmtPlaneVol->AddNode(pmtContVol, 20, trPmtCont4);
+ pmtPlaneVol->AddNode(pmtContVol, 23, trPmtCont5);
+ pmtPlaneVol->AddNode(pmtContVol, 26, trPmtCont6);
+
+ TGeoTranslation* trPmtPlane =
+ new TGeoTranslation(pmtXShift, 0., boxZLen / 2. - boxBackThickness - pmtZPadding + (pmtThickness / 2.));
+ gasVol->AddNode(pmtPlaneVol, 1, trPmtPlane);
+
+
+ gGeoMan->CheckOverlaps();
+ gGeoMan->PrintOverlaps();
+
+ //Draw
+ pmtPixelVol->SetLineColor(kYellow + 4);
+ aerogelVol->SetLineColor(kCyan);
+ aerogelVol->SetTransparency(70);
+ gasVol->SetLineColor(kGreen);
+ gasVol->SetTransparency(60);
+ kaptonVol->SetLineColor(kBlue);
+ kaptonVol->SetTransparency(70);
+
+
+ //gGeoMan->SetTopVisible();
+ //boxVol->SetVisibility(true);
+
+ // boxVol->Draw("ogl");
+ //top->Draw("ogl");
+ gGeoMan->SetVisLevel(7); //7
+
+ rich->Export(geoFileName); // an alternative way of writing the trd volume
+ TFile* geoFile = new TFile(geoFileName, "UPDATE");
+ trBox->Write();
+ cout << endl << "Geometry exported to " << geoFileName << endl;
+ geoFile->Close();
+}
diff --git a/macro/mcbm/mcbm_transport.C b/macro/mcbm/mcbm_transport.C
index 9963d13ad7..ca1be9118d 100644
--- a/macro/mcbm/mcbm_transport.C
+++ b/macro/mcbm/mcbm_transport.C
@@ -18,13 +18,14 @@
void SetTrack(CbmTransport*, Double_t, Int_t, Double_t, Double_t, Double_t);
-void mcbm_transport(Int_t nEvents = 10,
+void mcbm_transport(Int_t nEvents = 2,
// const char* setupName = "mcbm_beam_2022_06",
// const char* setupName = "mcbm_beam_2022_05",
// const char* setupName = "mcbm_beam_2022_04",
// const char* setupName = "mcbm_beam_2022_03",
// const char* setupName = "mcbm_beam_2022_02",
// const char* setupName = "mcbm_beam_2022_01",
+ // const char* setupName = "mcbm_beam_2022_03_09_carbon",
const char* setupName = "mcbm_beam_2021_07_surveyed",
// const char* setupName = "mcbm_beam_2021_04",
// const char* setupName = "mcbm_beam_2021_03",
@@ -32,7 +33,7 @@ void mcbm_transport(Int_t nEvents = 10,
// const char* setupName = "mcbm_beam_2019_11",
// const char* setupName = "mcbm_beam_2019_03",
// const char* setupName = "sis18_mcbm_25deg_long",
- const char* output = "data/test", const char* inputFile = "", Bool_t overwrite = kFALSE)
+ const char* output = "data/test", const char* inputFile = "", Bool_t overwrite = kTRUE)
{
// --- Define the beam angle ----------------------------------------------
Double_t beamRotY = 25.;
--
GitLab