diff --git a/macro/geometry/README.md b/macro/geometry/README.md
index 04114da6420b1ea3c3e1797579670bf23f584e8c..dc4f35f75c140ebc3de23bb3d055e537bc45e1e0 100644
--- a/macro/geometry/README.md
+++ b/macro/geometry/README.md
@@ -203,3 +203,21 @@ Allows to run the `check_radlen` script on sets of geometries present in the sub
Example for a version: `./macro/geometry/check_radlen_bulk.sh v22a` \
Special tags which cannot be used in this way: `git`, `all`, `main`
+
+---
+
+## examine_materials.C
+
+Root macro which compares the materials in the simulated TGeoManager to the materials as extracted from the geometry
+binaries for a specified setup. Compares only the names.
+
+Sugggested commands
+
+root -q examine_materials.C
+
+root -q 'examine_materials.C("sis100_muon_jpsi_DEC21")'
+
+ for file in `ls ../../geometry/setup/setup_mcbm_beam_2022_*`; do tag=`echo $file | awk -F'/' '{print $5}' | sed 's|setup_|"|' | sed 's_.C$_"_'`; echo $tag; root -q './examine_material.C('$tag')'; done > MCBM_2022_EXAMINE
+
+
+
diff --git a/macro/geometry/examine_material.C b/macro/geometry/examine_material.C
new file mode 100644
index 0000000000000000000000000000000000000000..4abc05100752efcc64f9b1b825e1939e8922c4b6
--- /dev/null
+++ b/macro/geometry/examine_material.C
@@ -0,0 +1,267 @@
+/* Copyright (C) 2023 Facility for Antiproton and Ion Research in Europe, Darmstadt
+ SPDX-License-Identifier: GPL-3.0-only
+ Authors: Eoin Clerkin [committer] */
+
+#define MAX_MATERIALS 200
+
+#include <cstdlib>
+#include <fstream>
+#include <iostream>
+
+#include <libgen.h>
+
+#define NODE_DEPTH 25
+
+/* The materials in the transported TGeoManager will be compared to the materials in the subsystem binaries.
+ The general idea pre-writing will be to run transport for a user specified setup generating the TGeoManager.
+ The number of volumes for each material will be recorded in a standard map. Each subsystem will be extracted
+ from the CBMSETUP, and the materials will be extraced from the respective binary and then substracted from
+ the from the TGeoManager version original. After going through all subsystems, if there is no volumes in the
+ original standard, the two TGeoManager and the material in the subsystem match, otherwise there is a material
+ mismatch detected. Output to the screen will help diagonsis. */
+
+char* padd(char* source_string, unsigned int size)
+{
+
+ char* padded_string = (char*) malloc(sizeof(char) * (size + 1));
+ strncpy(padded_string, source_string, size + 1);
+
+ unsigned int i = strlen(source_string);
+ if (i < size)
+ for (; i < size; i++)
+ padded_string[i] = ' ';
+ padded_string[size] = '\0';
+
+ return padded_string;
+};
+
+typedef struct MAT_VOL mat_vol_t;
+
+struct MAT_VOL {
+ int count; // Number of volumes
+ char name[100];
+};
+
+int REG_MATERIALS = 0; // counter for the currently number of registered materials.
+
+int get_register_count_material(const char* MATERIAL_NAME, mat_vol_t* material_array)
+{
+ int i = 0;
+
+ while (i <= MAX_MATERIALS && i <= REG_MATERIALS) {
+ // Register
+ if (REG_MATERIALS == i) {
+ strcpy((material_array[i]).name, MATERIAL_NAME);
+ REG_MATERIALS++;
+ material_array[i].count = 0;
+ break;
+ };
+ // Found
+ if (!strcmp((material_array[i]).name, MATERIAL_NAME)) { break; };
+ i++;
+ };
+ material_array[i].count++;
+
+ return 0;
+}
+
+int extract_mat_vol(char* fileName, mat_vol_t* material_array)
+{
+ TGeoNode* node[NODE_DEPTH];
+ TGeoMedium* med;
+ TGeoManager* gman = nullptr;
+ TGeoVolume* top = nullptr;
+
+ TFile file(fileName, "READ");
+ if (!file.IsOpen()) {
+ std::cerr << "Was a root geo file supplied?" << endl;
+ exit(1);
+ };
+
+ file.GetListOfKeys()->Print();
+ file.GetSize();
+ char geo_tag[100];
+ strcpy(geo_tag, fileName);
+ char* ptr = strstr(geo_tag, ".geo.root");
+ *(ptr) = '\0';
+
+ gman = (TGeoManager*) file.Get("FAIRGeom");
+ if (gman != nullptr) { top = gman->GetTopVolume(); };
+ if (top == nullptr) top = (TGeoVolume*) file.Get(geo_tag);
+ if (top == nullptr) top = (TGeoVolume*) file.Get("top");
+ if (top == nullptr) top = (TGeoVolume*) file.Get("TOP");
+ if (top == nullptr) top = (TGeoVolume*) file.Get("Top");
+ if (top == nullptr) top = (TGeoVolume*) file.Get("geometryFromGDML");
+ if (top == nullptr) top = (TGeoVolume*) file.Get(file.GetListOfKeys()->First()->GetName());
+ if (top == nullptr) {
+ std::cerr << "No Top Volume found. Is the TGeoManager or Top Volume unusally named" << std::endl;
+ exit(1);
+ };
+
+ TObjArray* nodes = top->GetNodes();
+
+ int i_array[NODE_DEPTH], num_array[NODE_DEPTH], j;
+ for (int i = 0; i < NODE_DEPTH; i++)
+ i_array[i] = 0;
+ for (int i = 0; i < NODE_DEPTH; i++)
+ num_array[i] = 0;
+
+ j = 0;
+ i_array[0] = 0;
+ num_array[0] = nodes->GetEntries();
+
+ while (num_array[0] != 0) {
+ if (j == 0) { node[0] = static_cast<TGeoNode*>(nodes->At(i_array[0])); }
+ else {
+ node[j] = static_cast<TGeoNode*>(node[j - 1]->GetDaughter(i_array[j]));
+ };
+ i_array[j]++; // Update number.
+ med = node[j]->GetMedium();
+ num_array[j + 1] = node[j]->GetNdaughters();
+
+ get_register_count_material(med->GetMaterial()->GetName(), material_array);
+
+ if (num_array[j + 1] > 0) {
+ j++;
+ num_array[j + 2] = 0;
+ };
+
+ while (i_array[j] == num_array[j]) {
+ num_array[j] = 0;
+ i_array[j] = 0; // Reset the counter before falling back.
+ j--;
+ };
+
+ if (j >= NODE_DEPTH) {
+ std::cerr << "Maximum depth of " << NODE_DEPTH << " exceeded. Increase NODE_DEPTH in header of macro," << endl;
+ exit(NODE_DEPTH);
+ };
+ };
+
+ file.Close();
+ return 0;
+}
+
+int examine_material(const char* setup = "sis100_electron")
+{
+ TString srcDir = gSystem->Getenv("VMCWORKDIR");
+ ECbmEngine engine = kGeant3;
+
+ CbmTransport run;
+ run.SetEngine(engine);
+ run.AddInput(srcDir + "/input/urqmd.auau.10gev.centr.root", kUnigen);
+ run.SetOutFileName("examine.tra.root", kTRUE);
+ run.SetParFileName("examaine.par.root");
+ run.SetGeoFileName("examine.geo.root");
+ run.LoadSetup(setup);
+ run.SetBeamPosition(0, 0, 0.1, 0.1);
+ run.SetRandomSeed(1234); // I fix this number so everything is deterministic
+ run.Run(1);
+
+ mat_vol_t* transport_materials = (mat_vol_t*) malloc(sizeof(mat_vol_t) * MAX_MATERIALS);
+ for (int i = 0; i < MAX_MATERIALS; i++) {
+ transport_materials[i].count = 0;
+ strcpy(transport_materials[i].name, "");
+ };
+
+ extract_mat_vol("examine.geo.root", transport_materials);
+
+ // Detector and passive subsystems
+ for (int i = 0; i < 120; i++)
+ printf("#");
+ printf("\n");
+
+ TString fileName;
+ char fileNameC[100];
+ CbmSetup* cs = run.GetSetup();
+ mat_vol_t* detector_materials = (mat_vol_t*) malloc(sizeof(mat_vol_t) * 500);
+ for (int i = 0; i < MAX_MATERIALS; i++) {
+ detector_materials[i].count = 0;
+ strcpy(detector_materials[i].name, "");
+ };
+
+ REG_MATERIALS = 0;
+ for (ECbmModuleId subsystem = ECbmModuleId::kMvd; subsystem < ECbmModuleId::kCave; ++subsystem) {
+ if (cs->IsActive(subsystem)) {
+
+ cs->GetGeoFileName(subsystem, fileName); // returns true
+
+
+ // Handling double beampipe
+ if ((subsystem == ECbmModuleId::kPipe) && strchr(fileName.Data(), ':')) {
+
+ fileName = srcDir + "/geometry/" + fileName; // Get the right directory
+ printf("fileName is %s\n", fileName.Data());
+ strcpy(fileNameC, fileName.Data());
+ char* change = strchr(fileNameC, ':');
+ *change = '\0';
+ printf("FIRST PIPE fileName now is %s\n", fileNameC);
+ extract_mat_vol(fileNameC, detector_materials);
+ strcpy(fileNameC, srcDir + "/geometry/" + (change + 1));
+ printf("SECOND PIPE fileName now is %s\n", fileNameC);
+ extract_mat_vol(fileNameC, detector_materials);
+
+ // Other subystems
+ }
+ else {
+
+ fileName = srcDir + "/geometry/" + fileName; // Get the right directory
+ printf("fileName is %s\n", fileName.Data());
+ strcpy(fileNameC, fileName.Data());
+ extract_mat_vol(fileNameC, detector_materials);
+ };
+ };
+ }
+
+ // Quick Sort Technique
+ mat_vol_t temp_mat;
+ for (int i = 0; i < (MAX_MATERIALS - 1); ++i) {
+ for (int j = 0; j < (MAX_MATERIALS - 1 - i); ++j) {
+ if ((detector_materials[j].count) < (detector_materials[j + 1].count)) {
+ temp_mat = detector_materials[j + 1];
+ detector_materials[j + 1] = detector_materials[j];
+ detector_materials[j] = temp_mat;
+ };
+ };
+ };
+ for (int i = 0; i < (MAX_MATERIALS - 1); ++i) {
+ for (int j = 0; j < (MAX_MATERIALS - 1 - i); ++j) {
+ if ((transport_materials[j].count) < (transport_materials[j + 1].count)) {
+ temp_mat = transport_materials[j + 1];
+ transport_materials[j + 1] = transport_materials[j];
+ transport_materials[j] = temp_mat;
+ };
+ };
+ };
+ int SUMT = 0, SUMD = 0;
+ for (int j = 0; j < MAX_MATERIALS; j++) {
+ if ((transport_materials + j)->count > 0) SUMT += (transport_materials + j)->count;
+ if ((detector_materials + j)->count > 0) SUMD += (detector_materials + j)->count;
+ };
+ (void) printf("TOTAL NUMBER OF VOLUMES\tTRANSPORTED: %d\t DETECTORS: %d\n", SUMT, SUMD);
+
+ printf("%s\t\tTransported\t\tSubsystems\t\tDifference \n", padd("Material", 25));
+ int FAIL = 0; // Success
+ int DIFF = 0;
+
+ for (int i = 0; i < MAX_MATERIALS; i++) {
+ if (0 < (transport_materials + i)->count) {
+ printf("%s\t\t%d", padd((transport_materials + i)->name, 25), (transport_materials + i)->count);
+ for (int j = 0; j < MAX_MATERIALS; j++) {
+ if (!strcmp((transport_materials[i].name), (detector_materials[j].name)))
+ printf("\t\t%d", (detector_materials + i)->count);
+ };
+ DIFF = ((transport_materials + i)->count) - ((detector_materials + i)->count);
+ printf("\t\t%d\n", DIFF);
+ };
+ if ((strcmp(transport_materials[i].name, "dummy") != 0) && (strcmp(transport_materials[i].name, "air") != 0))
+ FAIL += abs(DIFF);
+ };
+
+ if (0 == FAIL) { (void) printf("[SUCCESS] Excluding air and dummy materials, no material errors were detected. \n"); }
+ else {
+ (void) printf("[FAILED] Excluding air and dummy materials, %d issues with materials were detected. \n", FAIL);
+ };
+
+ return FAIL; // Success (when FAIL=0)
+}