Skip to content
Snippets Groups Projects

Add initial geometry version (MUST v24a) of PAnda STrAw modules in mCBM for Feb 2025

Merged Add initial geometry version (MUST v24a) of PAnda STrAw modules in mCBM for Feb 2025
d.emschermann/cbmroot_geometry:mcbm_pasta_v24a into master
All threads resolved!
Merged David Emschermann requested to merge d.emschermann/cbmroot_geometry:mcbm_pasta_v24a into master
All threads resolved!
Compare
and
3 files
+ 408
− 0
Compare changes
  • Side-by-side
  • Inline
Files
3
macro/must/mcbm/Create_MUST_Geometry_v24a.C 0 → 100644
+ 393
− 0
/* Copyright (C) 2020-2025 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt
SPDX-License-Identifier: GPL-3.0-only
Authors: Florian Uhlig, David Emschermann [committer] */
// 2024-12-12 - v24a - DE - first version for PAnda STrAws in mCBM, 4 modules
// 2024-11-25 - v01a - RK - first version
// details concering the copynr of layers, modules and straw assemblies
// 1 0 4 0 1 0 6 4 sample copy number of a straw
// | | / | / \---- number of straw in module (1-128, 3 decimal digits, 7 bits -127)
// | | \-------- number of module in layer (1-36, 2 decimal digits, 6 bits - 63)
// | \------------ number of layer in MuST (1-12, 2 decimal digits, 4 bits - 15)
// \-------------- leading 1 to allow for 0s in layer number
// clang-format off
#include <iostream>
#include "TGeoBBox.h"
#include "TGeoManager.h"
#include "TGeoVolume.h"
void dump_info_file();
using namespace std;
// Name of output file with geometry
const TString tagVersion = "v24a";
// const TString geoVersion = "must_" + tagVersion;
const TString geoVersion = "must_" + tagVersion + "_mcbm";
const TString FileNameSim = geoVersion + ".geo.root?reproducible"; // Reproducible flag removes changing timestamps etc...
const TString FileNameGeo = geoVersion + "_geo.root";
const TString FileNameInfo = geoVersion + ".geo.info";
// Double_t fTranslationZ = 0.; // cm - Z Translation of the whole Geometry
Double_t fTranslation[3] = { 50., 0., 340. + 1.5 }; // cm - Translation of the whole Geometry
// DE-1layer const Int_t nLayer = 1; // number of layers
// DE-1layer Double_t StereoAngle[nLayer] = { 0. }; // X,U,V,X stereo angle of straws
// DE-1layer const Int_t mustlayout[nLayer][nModule] = { { 1 } };
// mCBM 2025
const Int_t nLayer = 4; // number of layers
const Int_t nModule = 1; // number of modules
const Int_t nModuleType = 4; // number of modules types
const Double_t zLayerPitch = 4.0; // distance of module front to module front in z
Double_t StereoAngle[nLayer] = { 0., -5., 5., 0.}; // X,U,V,X stereo angle of straws - mCBM 2025
// Double_t StereoAngle[nLayer] = { 0., 0., 0., 0.}; // X,U,V,X stereo angle of straws
// Double_t StereoAngle[nLayer] = { 0., 5., -5., 0.}; // X,U,V,X stereo angle of straws
const Int_t mustlayout[nLayer][nModule] = { { 1 }, { 1 }, { 1 }, { 1 } };
// mCBM ------------
// DE4 // CBM FAIR
// DE4 const Int_t nLayer = 12; // number of layers
// DE4 const Int_t nModule = 36; // number of modules
// DE4
// DE4 const Int_t mustlayout[nLayer][nModule] =
// DE4 // station 1
// DE4 { { 1, 1, 1, 1, 1, 1, 1, 2, 4, 3, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 4, 3, 2, 1, 1, 1, 1, 1, 1, 1 },
// DE4 { 1, 1, 1, 1, 1, 1, 1, 2, 4, 3, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 4, 3, 2, 1, 1, 1, 1, 1, 1, 1 },
// DE4 { 1, 1, 1, 1, 1, 1, 1, 2, 4, 3, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 4, 3, 2, 1, 1, 1, 1, 1, 1, 1 },
// DE4 { 1, 1, 1, 1, 1, 1, 1, 2, 4, 3, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 4, 3, 2, 1, 1, 1, 1, 1, 1, 1 },
// DE4 // station 2
// DE4 { 1, 1, 1, 1, 1, 1, 1, 2, 4, 3, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 4, 3, 2, 1, 1, 1, 1, 1, 1, 1 },
// DE4 { 1, 1, 1, 1, 1, 1, 1, 2, 4, 3, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 4, 3, 2, 1, 1, 1, 1, 1, 1, 1 },
// DE4 { 1, 1, 1, 1, 1, 1, 1, 2, 4, 3, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 4, 3, 2, 1, 1, 1, 1, 1, 1, 1 },
// DE4 { 1, 1, 1, 1, 1, 1, 1, 2, 4, 3, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 4, 3, 2, 1, 1, 1, 1, 1, 1, 1 },
// DE4 // station 3
// DE4 { 1, 1, 1, 1, 1, 1, 1, 2, 4, 3, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 4, 3, 2, 1, 1, 1, 1, 1, 1, 1 },
// DE4 { 1, 1, 1, 1, 1, 1, 1, 2, 4, 3, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 4, 3, 2, 1, 1, 1, 1, 1, 1, 1 },
// DE4 { 1, 1, 1, 1, 1, 1, 1, 2, 4, 3, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 4, 3, 2, 1, 1, 1, 1, 1, 1, 1 },
// DE4 { 1, 1, 1, 1, 1, 1, 1, 2, 4, 3, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 4, 3, 2, 1, 1, 1, 1, 1, 1, 1 } };
// DE5 // conical beampipe support in S2 and S3 modules
// DE5 const Int_t mustlayout[nLayer][nModule] =
// DE5 // station 1
// DE5 { { 1, 1, 1, 1, 1, 1, 1, 2, 4, 3, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 4, 3, 2, 1, 1, 1, 1, 1, 1, 1 },
// DE5 { 1, 1, 1, 1, 1, 1, 1, 2, 4, 3, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 4, 3, 2, 1, 1, 1, 1, 1, 1, 1 },
// DE5 { 1, 1, 1, 1, 1, 1, 1, 2, 4, 3, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 4, 3, 2, 1, 1, 1, 1, 1, 1, 1 },
// DE5 { 1, 1, 1, 1, 1, 1, 1, 2, 4, 3, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 4, 3, 2, 1, 1, 1, 1, 1, 1, 1 },
// DE5 // station 2
// DE5 { 1, 1, 1, 1, 1, 1, 1, 2, 5, 6, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 5, 6, 2, 1, 1, 1, 1, 1, 1, 1 },
// DE5 { 1, 1, 1, 1, 1, 1, 1, 2, 5, 6, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 5, 6, 2, 1, 1, 1, 1, 1, 1, 1 },
// DE5 { 1, 1, 1, 1, 1, 1, 1, 2, 5, 6, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 5, 6, 2, 1, 1, 1, 1, 1, 1, 1 },
// DE5 { 1, 1, 1, 1, 1, 1, 1, 2, 5, 6, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 5, 6, 2, 1, 1, 1, 1, 1, 1, 1 },
// DE5 // station 3
// DE5 { 1, 1, 1, 1, 1, 1, 1, 2, 7, 8, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 7, 8, 2, 1, 1, 1, 1, 1, 1, 1 },
// DE5 { 1, 1, 1, 1, 1, 1, 1, 2, 7, 8, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 7, 8, 2, 1, 1, 1, 1, 1, 1, 1 },
// DE5 { 1, 1, 1, 1, 1, 1, 1, 2, 7, 8, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 7, 8, 2, 1, 1, 1, 1, 1, 1, 1 },
// DE5 { 1, 1, 1, 1, 1, 1, 1, 2, 7, 8, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 7, 8, 2, 1, 1, 1, 1, 1, 1, 1 } };
// some global variables
TGeoManager* gGeoMan = NULL; // Pointer to TGeoManager instance
TGeoVolume* gModules[nModuleType]; // Global storage for module types
// Forward declarations
void create_materials_from_media_file();
TGeoVolume* create_must_module_type(Int_t moduleType, Int_t modulecopynr);
void Create_MUST_Geometry_v24a()
{
// Load FairRunSim to ensure the correct unit system
FairRunSim* sim = new FairRunSim();
// Load needed material definition from media.geo file
create_materials_from_media_file();
// Get the GeoManager for later usage
gGeoMan = (TGeoManager*) gROOT->FindObject("FAIRGeom");
gGeoMan->SetVisLevel(10);
// Create the top volume
// TGeoBBox* topbox = new TGeoBBox("", 1000., 1000., 2000.);
TGeoBBox* topbox = new TGeoBBox("", 100., 1000., 2000.);
TGeoVolume* top = new TGeoVolume("top", topbox, gGeoMan->GetMedium("air"));
gGeoMan->SetTopVolume(top);
TGeoVolume* must = new TGeoVolumeAssembly(geoVersion);
top->AddNode(must, 1);
for (Int_t iLayer=1; iLayer <= nLayer; iLayer++)
{
// add layer keeping volume
TString layername = Form("layer%02d", iLayer);
TGeoVolume* layer = new TGeoVolumeAssembly(layername);
Int_t layercopynr = 100 + iLayer;
// cout << "layercopynr " << layercopynr << endl;
must->AddNode(layer, layercopynr);
for (Int_t iModule=1; iModule <= nModule; iModule++)
{
Int_t modulecopynr = layercopynr * 100 + iModule;
// cout << "modulecopynr " << modulecopynr << endl;
gModules[iModule] = create_must_module_type(iModule, modulecopynr);
TGeoRotation* module_rotation = new TGeoRotation();
module_rotation->RotateZ(StereoAngle[iLayer-1]);
TGeoCombiTrans* module_placement = new TGeoCombiTrans(0., 0., (iLayer-1) * zLayerPitch, module_rotation);
layer->AddNode(gModules[iModule], modulecopynr, module_placement);
}
}
gGeoMan->CloseGeometry();
gGeoMan->CheckOverlaps(0.001);
gGeoMan->PrintOverlaps();
gGeoMan->Test();
must->Export(FileNameSim); // an alternative way of writing the must volume
TFile* outfile = new TFile(FileNameSim, "UPDATE");
TGeoTranslation* must_placement =
new TGeoTranslation("must_trans", fTranslation[0], fTranslation[1], fTranslation[2]);
must_placement->Write();
outfile->Close();
outfile = new TFile(FileNameGeo, "RECREATE");
gGeoMan->Write(); // use this is you want GeoManager format in the output
outfile->Close();
dump_info_file();
top->Draw("ogl");
}
void create_materials_from_media_file()
{
// Use the FairRoot geometry interface to load the media which are already defined
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();
// Read the required media and create them in the GeoManager
FairGeoMedia* geoMedia = geoFace->getMedia();
FairGeoBuilder* geoBuild = geoLoad->getGeoBuilder();
//Media
FairGeoMedium *mcarbon = geoMedia->getMedium("TRDcarbon");
if (mcarbon == NULL) Fatal("Main", "FairMedium TRDcarbon not found");
geoBuild->createMedium(mcarbon);
FairGeoMedium *mkapton = geoMedia->getMedium("MUCHkapton");
if (mkapton == NULL) Fatal("Main", "FairMedium MUCHkapton not found");
geoBuild->createMedium(mkapton);
FairGeoMedium *mMUCHargon_co2 = geoMedia->getMedium("MUCHGEMmixture");
if (mMUCHargon_co2 == NULL) Fatal("Main", "FairMedium mMUCHargon_co2 not found");
geoBuild->createMedium(mMUCHargon_co2);
//---------------------------------------------------------------------------------
// FairGeoMedium* mAluminium = geoMedia->getMedium("aluminium");
// if (mAluminium == NULL) Fatal("Main", "FairMedium aluminium not found");
// geoBuild->createMedium(mAluminium);
//
// FairGeoMedium* mMylar = geoMedia->getMedium("mylar");
// if (mMylar == NULL) Fatal("Main", "FairMedium mylar not found");
// geoBuild->createMedium(mMylar);
//
// FairGeoMedium* mIron = geoMedia->getMedium("iron"); // TODO should be tungsten?
// if (mIron == NULL) Fatal("Main", "FairMedium iron not found");
// geoBuild->createMedium(mIron);
//
// FairGeoMedium* mAir = geoMedia->getMedium("air");
// if (mAir == NULL) Fatal("Main", "FairMedium air not found");
// geoBuild->createMedium(mAir);
}
TGeoVolume* create_must_module_type(Int_t moduleType, Int_t modulecopynr)
{
TString modulename = Form("module%02d", moduleType);
TGeoVolume* module = new TGeoVolumeAssembly(modulename);
//------------------------------------------------------------------
TGeoMedium* carbon = gGeoMan->GetMedium("TRDcarbon"); // Carbon
TGeoMedium* kapton = gGeoMan->GetMedium("MUCHkapton"); // Kapton
TGeoMedium* argon_co2 = gGeoMan->GetMedium("MUCHGEMmixture"); // active medium (Ar + CO2)
//------------------------------------------------------------------
// not used TGeoMedium* medAl = gGeoMan->GetMedium("aluminium");
// not used if (medAl == NULL) Fatal("Main", "Medium vacuum not found");
// not used
// not used TGeoMedium* medMylar = gGeoMan->GetMedium("mylar");
// not used if (medMylar == NULL) Fatal("Main", "Medium mylar not found");
// not used
// not used TGeoMedium* medIron = gGeoMan->GetMedium("iron");
// not used if (medIron == NULL) Fatal("Main", "Medium iron not found");
// not used
// not used TGeoMedium* medAir = gGeoMan->GetMedium("air");
// not used if (medAir == NULL) Fatal("Main", "Medium air not found");
// Straw tube dimensions
Double_t strawLength = 230.0; // cm
Double_t strawRadius = 0.245; // cm
Double_t wallThickness = 0.0065; // cm
// Module and straw layout parameters
Int_t nStrawsPerColumnLong = 64; // DE
Int_t nStrawsPerModuleS3 = 32; // DE
Double_t spacingY = 0.550; // cm // DE
Double_t spacingX = 0.525; // cm // DE
Double_t offsetX = spacingX / 2.0;
//Double_t size_module_long = std::ceil(spacingX * nStrawsPerColumnLong); // DE
Double_t xwidth_module_long = 34.0; // cm // DE
// (34 / 2. - 31.75 * 0.525 + 15.75 * 0.525) * 2 = 17.2
Double_t xwidth_module_S3 = 17.2; // cm // DE
Double_t moduleXSize = xwidth_module_S3;
Double_t moduleYSize = strawLength;
Double_t moduleZSize = 3.;
Double_t frameThick = 0.05;
Double_t moduleXPos = 0., moduleYPos = 0., moduleZPos = moduleZSize/2.;
// Carbon SideFrame
TGeoBBox* must_vert_frame = new TGeoBBox("must_vert_frame", frameThick/2., moduleYSize/2., moduleZSize/2.); // vertical
TGeoVolume* must_vert_frame_vol = new TGeoVolume("framevert", must_vert_frame, carbon);
must_vert_frame_vol->SetLineColor(kBlue);
TGeoTranslation* frameTransPosX = new TGeoTranslation(moduleXPos + (moduleXSize/2. - frameThick/2.), moduleYPos, moduleZPos);
module->AddNode(must_vert_frame_vol, 1, frameTransPosX);
TGeoTranslation* frameTransNegX = new TGeoTranslation(moduleXPos - (moduleXSize/2. - frameThick/2.), moduleYPos, moduleZPos);
module->AddNode(must_vert_frame_vol, 2, frameTransNegX);
TGeoVolumeAssembly* strawAssembly = new TGeoVolumeAssembly("strawAssembly");
TGeoRotation* strawRotat = new TGeoRotation("strawRotat", 0.0, 90.0, 0.0);
// Common straw geometry
TGeoTube* strawGas = new TGeoTube("StrawGas", 0, strawRadius - wallThickness, strawLength / 2.);
TGeoVolume* volStrawGas = new TGeoVolume("strawgasactive", strawGas, argon_co2);
volStrawGas->SetLineColor(kGreen);
volStrawGas->SetTransparency(40);
strawAssembly->AddNode(volStrawGas,1);
TGeoTube* strawWall = new TGeoTube("StrawWall", strawRadius - wallThickness, strawRadius, strawLength / 2.);
TGeoVolume* volStrawWall = new TGeoVolume("strawwallpassive", strawWall, kapton);
volStrawWall->SetLineColor(kRed);
//volStrawWall->SetTransparency(40);
strawAssembly->AddNode(volStrawWall,1);
// not used double strawRad = 0.25, straw_Length = moduleYSize;
// not used double stGasRad = 0.24, stGas_Length = straw_Length;
// not used double wireRad = 0.01, wire_Length = straw_Length;
// not used
// not used auto wireShape = new TGeoTube("wireShape", 0., wireRad, wire_Length/2.);
// not used TGeoVolume* wireVol = new TGeoVolume("MustWire", wireShape, medIron);
// not used wireVol->SetLineColor(kRed);
// not used strawAssembly->AddNode(wireVol,1);
// not used
// not used auto stGasShape = new TGeoTube("stGasShape", wireRad, stGasRad, stGas_Length/2.);
// not used TGeoVolume* stGasVol = new TGeoVolume("MustStGas", stGasShape, medAir);
// not used stGasVol->SetTransparency();
// not used stGasVol->SetLineColor(kYellow);
// not used strawAssembly->AddNode(stGasVol,1);
// not used
// not used auto strawShape = new TGeoTube("strawShape", stGasRad, strawRad, straw_Length/2.);
// not used TGeoVolume* strawVol = new TGeoVolume("MustStraw", strawShape, medMylar);
// not used strawVol->SetLineColor(kGreen);
// not used strawAssembly->AddNode(strawVol,1);
Int_t strawcopynr = 0;
Int_t nStraws = 2 * nStrawsPerModuleS3;
for ( int istraw = 0; istraw < nStraws; istraw++ )
{
Double_t xoffset = 15.75;
Double_t yoffset = -spacingY/2.;
if (istraw%2==1)
{
xoffset = 15.25;
yoffset = spacingY/2.;
}
TGeoTranslation* strawTrans = new TGeoTranslation((istraw/2 - xoffset) * spacingX, 0., yoffset + moduleZPos);
TGeoCombiTrans* strawCombi = new TGeoCombiTrans(*strawTrans, *strawRotat);
strawcopynr = modulecopynr * 1000 + istraw + 1;
// cout << "strawcopynr " << strawcopynr << endl;
module->AddNode(strawAssembly, strawcopynr, strawCombi);
// TGeoTranslation* strawTrans = new TGeoTranslation((istraw/2 - 15.75) * spacingX, 0., -spacingY/2.);
// TGeoTranslation* strawTrans = new TGeoTranslation((istraw/2 - 15.25) * spacingX, 0., spacingY/2.);
}
// for ( int istraw = 0 ; istraw < nStrawsPerModuleS3 ; istraw++ )
// {
// TGeoTranslation* strawTrans = new TGeoTranslation((istraw - 15.75) * spacingX, 0., -spacingY/2.);
// TGeoCombiTrans* strawCombi = new TGeoCombiTrans(*strawTrans, *strawRotat);
//
// strawcopynr = modulecopynr * 1000 + istraw + 1;
// // cout << "strawcopynr " << strawcopynr << endl;
// module->AddNode(strawAssembly, strawcopynr, strawCombi);
// }
//
// for ( int istraw = 0 ; istraw < nStrawsPerModuleS3 ; istraw++ ) {
// TGeoTranslation* strawTrans = new TGeoTranslation((istraw - 15.25) * spacingX, 0., spacingY/2.);
// TGeoCombiTrans* strawCombi = new TGeoCombiTrans(*strawTrans, *strawRotat);
//
// strawcopynr = modulecopynr * 1000 + istraw + 1 + 32;
// // cout << "strawcopynr " << strawcopynr << endl;
// module->AddNode(strawAssembly, strawcopynr, strawCombi);
// }
return module;
}
void dump_info_file()
{
TDatime datetime; // used to get timestamp
printf("writing summary information file: %s\n", FileNameInfo.Data());
FILE* ifile;
ifile = fopen(FileNameInfo.Data(), "w");
if (ifile == NULL) {
printf("error opening %s\n", FileNameInfo.Data());
exit(1);
}
fprintf(ifile, "#\n## %s information file\n#\n\n", geoVersion.Data());
fprintf(ifile, "# created %d\n\n", datetime.GetDate());
fprintf(ifile, "# MUST geometry consisting fo %d layers\n", nLayer);
for (Int_t iLayer=1; iLayer <= nLayer; iLayer++)
{
fprintf(ifile, "Layer %d: number of modules: %d module type %d StereoAngle %f.1 degrees\n", iLayer, nModule, nModuleType, StereoAngle[iLayer]);
}
fprintf(ifile, "Front module x-offset %f.1 cm\n", fTranslation[0]);
fprintf(ifile, "Front module y-offset %f.1 cm\n", fTranslation[1]);
fprintf(ifile, "Front module z-offset %f.1 cm\n", fTranslation[2]);
fprintf(ifile, "\n");
fclose(ifile);
}