// -------------------------------------------------------------------------- // // Macro for standard transport simulation with GEEANT3 // in the CBM muon setup: STS + MUCH + TRD + TOF // // No MVD nor ECAL // // Input is UrQMD // // Parameters are stored in the output file. // // A.Kiseleva // Add-ons, cosmetics and comments V.Friese 2007-11-09 // -------------------------------------------------------------------------- void much_sim(Int_t nEvents = 10000) { // ======================================================================== // Adjust this part according to your requirements // Int_t iFirst = 1; TString system = "auau"; TString beam = "25gev"; TString trigger = "centr"; TString particle = "omega"; //TString directory = "/d/cbm02/andrey/events/much/compact/signal/"; //TString dir = "/d/cbm02/andrey/events/much/10stations/"; TString dir = "/d/cbm02/andrey/events/much/10stations/signal/"; TString outFile = dir + beam + "/" + particle + "/mc." + system + "." + beam + "." + particle + "." + trigger + ".root"; TString parFile = dir + beam + "/" + particle + "/params." + system + "." + beam + "." + particle + "." + trigger + ".root"; TString inFile = "/d/cbm03/urqmd/" + system + "/" + beam + "/" + trigger + "/urqmd." + system + "." + beam + "." + trigger+ ".0001.ftn14"; TString plutoFile = "/u/andrey/cbm/much/pluto/" + particle + "/" + beam + "/" + particle + "." + "0000" + ".root"; // ----- Confirm input parameters ------------------------------------ cout << endl; cout << "======== CBMROOT Macro much_sim =================" << endl; cout << "Input file is " << inFile << endl; cout << "Output file is " << outFile << endl; cout << "Events to process: " << nEvents << endl; // cout << "First event is " << iFirst << endl; cout << "===================================================" << endl; cout << endl; // ----- Specify MUCH related geometry -------------------------------- // much_standard.geo: Full absorber system, 18 detectors // much_compact.geo: Without last absorber, 15 detectors // Use pipe_much.geo for the beam pipe in both cases. // In case you want the addtional W shielding around the pipe, // use shield_standard.geo or shield_compact.geo, respective to the // MUCH geometry. Otherwise, define an empty string. //TString muchGeom = "much_standard.geo"; // TString muchGeom = "../much/geometry/much_compact.geo"; TString muchGeom = "../much/geometry/much_compact_10_2.geo"; TString pipeGeom = "pipe_much.geo"; TString shieldGeom = "shield_standard.geo"; // ----- Other geometries --------------------------------------------- TString caveGeom = "cave.geo"; TString targetGeom = "target.geo"; TString magnetGeom = "magnet_muon.geo"; TString stsGeom = "sts_allstrips.geo"; TString trdGeom = "trd_standard.geo"; TString tofGeom = "tof_standard.geo"; // ----- Magnetic field ----------------------------------------------- TString fieldMap = "FieldMuonMagnet"; // name of field map Double_t fieldZ = 50.; // field centre z position Double_t fieldScale = 1.; // field scaling factor // In general, the following parts need not be touched // ======================================================================== // ----- Set unique random generator seed ----------------------------- // Comment this out if you want to have a defined seed for reproducability. gRandom->SetSeed(0); // ------------------------------------------------------------------------ // ---- Debug option ------------------------------------------------- gDebug = 0; // ------------------------------------------------------------------------ // ----- Timer -------------------------------------------------------- cout << endl << "=== much_sim.C : Starting timer ..." << endl; TStopwatch timer; timer.Start(); // ------------------------------------------------------------------------ // ---- Load libraries ------------------------------------------------- cout << endl << "=== much_sim.C : Loading libraries ..." << endl; gROOT->LoadMacro("$VMCWORKDIR/gconfig/basiclibs.C"); basiclibs(); gSystem->Load("libGeoBase"); gSystem->Load("libParBase"); gSystem->Load("libBase"); gSystem->Load("libCbmBase"); gSystem->Load("libField"); gSystem->Load("libGen"); gSystem->Load("libPassive"); gSystem->Load("libSts"); gSystem->Load("libMuch"); gSystem->Load("libTrd"); gSystem->Load("libTof"); // ----------------------------------------------------------------------- // ----- Create simulation run ---------------------------------------- cout << endl << "=== much_sim.C : Creating run and database ..." << endl; CbmRunSim* fRun = new CbmRunSim(); fRun->SetName("TGeant3"); // Transport engine fRun->SetOutputFile(outFile); // Output file CbmRuntimeDb* rtdb = fRun->GetRuntimeDb(); // ------------------------------------------------------------------------ // ----- Create media ------------------------------------------------- cout << endl << "=== much_sim.C : Set materials ..." << endl; fRun->SetMaterials("media.geo"); // Materials // ------------------------------------------------------------------------ // ----- Create detectors and passive volumes ------------------------- cout << endl << "=== much_sim.C : Create geeometry ..." << endl; if ( caveGeom != "" ) { CbmModule* cave = new CbmCave("CAVE"); cave->SetGeometryFileName(caveGeom); fRun->AddModule(cave); cout << " --- " << caveGeom << endl; } if ( pipeGeom != "" ) { CbmModule* pipe = new CbmPipe("PIPE"); pipe->SetGeometryFileName(pipeGeom); fRun->AddModule(pipe); cout << " --- " << pipeGeom << endl; } if ( shieldGeom != "" ) { CbmModule* shield = new CbmShield("SHIELD"); shield->SetGeometryFileName(shieldGeom); fRun->AddModule(shield); cout << " --- " << shieldGeom << endl; } if ( targetGeom != "" ) { CbmModule* target = new CbmTarget("Target"); target->SetGeometryFileName(targetGeom); fRun->AddModule(target); cout << " --- " << targetGeom << endl; } if ( magnetGeom != "" ) { CbmModule* magnet = new CbmMagnet("MAGNET"); magnet->SetGeometryFileName(magnetGeom); fRun->AddModule(magnet); cout << " --- " << magnetGeom << endl; } if ( stsGeom != "" ) { CbmDetector* sts = new CbmSts("STS", kTRUE); sts->SetGeometryFileName(stsGeom); fRun->AddModule(sts); cout << " --- " << stsGeom << endl; } if ( muchGeom != "" ) { CbmDetector* much = new CbmMuch("MUCH", kTRUE); much->SetGeometryFileName(muchGeom); fRun->AddModule(much); cout << " --- " << muchGeom << endl; } if ( trdGeom != "" ) { CbmDetector* trd = new CbmTrd("TRD",kTRUE ); trd->SetGeometryFileName(trdGeom); fRun->AddModule(trd); cout << " --- " << trdGeom << endl; } if ( tofGeom != "" ) { CbmDetector* tof = new CbmTof("TOF", kTRUE); tof->SetGeometryFileName(tofGeom); fRun->AddModule(tof); cout << " --- " << tofGeom << endl; } // ------------------------------------------------------------------------ // ----- Create magnetic field ---------------------------------------- cout << endl << "=== much_sim.C : Create magnetic field ..." << endl; CbmFieldMap* magField = NULL; if ( fieldMap == "FieldActive" || fieldMap == "FieldIron") magField = new CbmFieldMapSym3(fieldMap); else if ( fieldMap == "FieldAlligator" ) magField = new CbmFieldMapSym2(fieldMap); else if ( fieldMap = "FieldMuonMagnet" ) magField = new CbmFieldMapSym3(fieldMap); else { cout << "===> ERROR: Field map " << fieldMap << " unknown! " << endl; exit; } magField->SetPosition(0., 0., fieldZ); magField->SetScale(fieldScale); fRun->SetField(magField); cout << " --- Field is " << fieldMap << endl; cout << " --- Position is z = " << fieldZ << " cm" << endl; cout << " --- Scale is " << fieldScale << endl; // ------------------------------------------------------------------------ // ----- Create PrimaryGenerator -------------------------------------- cout << endl << "=== much_sim.C : Create generators ..." << endl; CbmPrimaryGenerator* primGen = new CbmPrimaryGenerator(); CbmPlutoGenerator *plutoGen= new CbmPlutoGenerator(plutoFile); primGen->AddGenerator(plutoGen); // CbmUrqmdGenerator* urqmdGen = new CbmUrqmdGenerator(inFile); // primGen->AddGenerator(urqmdGen); fRun->SetGenerator(primGen); // ------------------------------------------------------------------------ // ----- Run initialisation ------------------------------------------- cout << endl << "=== much_sim.C : Initialise run ..." << endl; fRun->Init(); // ------------------------------------------------------------------------ // ----- Runtime database --------------------------------------------- cout << endl << "=== much_sim.C : Set up database ..." << endl; cout << " Parameters will be saved to output file" << endl; CbmFieldPar* fieldPar = (CbmFieldPar*) rtdb->getContainer("CbmFieldPar"); fieldPar->SetParameters(magField); fieldPar->setChanged(); fieldPar->setInputVersion(fRun->GetRunId(),1); Bool_t kParameterMerged = kTRUE; CbmParRootFileIo* parOut = new CbmParRootFileIo(kParameterMerged); parOut->open(parFile.Data());//gFile); rtdb->setOutput(parOut); rtdb->saveOutput(); rtdb->print(); // ------------------------------------------------------------------------ // ----- Start run ---------------------------------------------------- cout << endl << "=== much_sim.C : Start run ..." << endl; fRun->Run(nEvents); // ------------------------------------------------------------------------ // ----- Finish ------------------------------------------------------- timer.Stop(); Double_t rtime = timer.RealTime(); Double_t ctime = timer.CpuTime(); cout << endl << endl; cout << "=== much_sim.C : Macro finished successfully." << endl; cout << "=== much_sim.C : Output file is " << outFile << endl; cout << "=== much_sim.C : Real time used: " << rtime << "s " << endl; cout << "=== much_sim.C : CPU time used : " << ctime << "s " << endl; cout << endl << endl; // ------------------------------------------------------------------------ }