MpdStrawECT.cxx 21 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532
  1. //------------------------------------------------------------------------------------------------------------------------
  2. // -------------------------------------------------------------------------
  3. // ----- MpdStrawECT source file -----
  4. // -------------------------------------------------------------------------
  5. #include <iostream>
  6. #include "TClonesArray.h"
  7. #include "TLorentzVector.h"
  8. #include "TParticle.h"
  9. #include "TVirtualMC.h"
  10. #include "TGeoManager.h"
  11. #include "TGeoMatrix.h"
  12. #include "FairGeoInterface.h"
  13. #include "FairGeoLoader.h"
  14. #include "FairGeoNode.h"
  15. #include "FairGeoRootBuilder.h"
  16. #include "MpdStack.h"
  17. #include "MpdStrawECTGeo.h"
  18. #include "FairRootManager.h"
  19. #include "MpdStrawECT.h"
  20. #include "MpdStrawECTPoint.h"
  21. #include "FairRuntimeDb.h"
  22. #include "MpdStrawECTGeoPar.h"
  23. #include "TObjArray.h"
  24. #include "FairRun.h"
  25. #include "FairVolume.h"
  26. #include "TMath.h"
  27. #include "TParticlePDG.h"
  28. #include "TGeoManager.h"
  29. #include "TGDMLParse.h"
  30. #include "FairGeoMedia.h"
  31. #include "TString.h"
  32. class FairVolume;
  33. //------------------------------------------------------------------------------------------------------------------------
  34. MpdStrawECT::MpdStrawECT()
  35. : FairDetector("strawECT", kTRUE)
  36. {
  37. fStrawECTCollection = new TClonesArray("MpdStrawECTPoint");
  38. fPosIndex = 0;
  39. fVerboseLevel = 1;
  40. }
  41. //------------------------------------------------------------------------------------------------------------------------
  42. MpdStrawECT::MpdStrawECT(const char* name, Bool_t active)
  43. : FairDetector(name, active)
  44. {
  45. fStrawECTCollection = new TClonesArray("MpdStrawECTPoint");
  46. fPosIndex = 0;
  47. fVerboseLevel = 1;
  48. }
  49. //------------------------------------------------------------------------------------------------------------------------
  50. MpdStrawECT::~MpdStrawECT()
  51. {
  52. if(fStrawECTCollection){ fStrawECTCollection->Delete(); delete fStrawECTCollection; }
  53. }
  54. //------------------------------------------------------------------------------------------------------------------------
  55. int MpdStrawECT::DistAndPoints(TVector3 p1, TVector3 p2, TVector3 p3, TVector3 p4, TVector3& pa, TVector3& pb) {
  56. TVector3 A = p2 - p1;
  57. TVector3 B = p4 - p3;
  58. TVector3 C = p1 - p3;
  59. if (p3 == p4) {
  60. pb = p4;
  61. TVector3 unit = A.Unit();
  62. Double_t dist = unit.Dot(-C);
  63. pa = p1 + dist * unit;
  64. return 0;
  65. }
  66. Double_t numer = C.Dot(B) * B.Dot(A) - C.Dot(A) * B.Dot(B);
  67. Double_t denom = A.Dot(A) * B.Dot(B) - B.Dot(A) * B.Dot(A);
  68. if (denom == 0) {
  69. // parallel lines
  70. pa.SetXYZ(-10000,-10000,-10000);
  71. pb.SetXYZ(-10000,-10000,-10000);
  72. return 1;
  73. }
  74. Double_t ma = numer / denom;
  75. Double_t mb = ( C.Dot(B) + B.Dot(A)*ma ) / B.Dot(B);
  76. pa = p1 + A*ma;
  77. pb = p3 + B*mb;
  78. return 0;
  79. }
  80. //------------------------------------------------------------------------------------------------------------------------
  81. TVector3 MpdStrawECT::GlobalToLocal(TVector3& global) {
  82. Double_t globPos[3];
  83. Double_t localPos[3];
  84. global.GetXYZ(globPos);
  85. gMC->Gmtod(globPos, localPos, 1);
  86. return TVector3(localPos);
  87. }
  88. //------------------------------------------------------------------------------------------------------------------------
  89. TVector3 MpdStrawECT::LocalToGlobal(TVector3& local) {
  90. Double_t globPos[3];
  91. Double_t localPos[3];
  92. local.GetXYZ(localPos);
  93. gMC->Gdtom(localPos, globPos, 1);
  94. return TVector3(globPos);
  95. }
  96. //----------------------------------------------------------------------------------------------------------------------
  97. Bool_t MpdStrawECT::ProcessHits(FairVolume* vol) {
  98. // Set parameters at entrance of volume. Reset ELoss.
  99. if(gMC->IsTrackEntering()) {
  100. ResetParameters();
  101. fELoss = 0.;
  102. fTime = gMC->TrackTime() * 1.0e09;
  103. fLength = gMC->TrackLength();
  104. fIsPrimary = 0;
  105. fCharge = -1;
  106. fPdgId = 0;
  107. TLorentzVector PosIn;
  108. gMC->TrackPosition(PosIn);
  109. fPosIn.SetXYZ(PosIn.X(), PosIn.Y(), PosIn.Z());
  110. gMC->TrackMomentum(fMom);
  111. TParticle* part = 0;
  112. part = gMC->GetStack()->GetCurrentTrack();
  113. if (part) {
  114. fIsPrimary = (Int_t)part->IsPrimary();
  115. fCharge = (Int_t)part->GetPDG()->Charge();
  116. fPdgId = (Int_t)part->GetPdgCode();
  117. }
  118. fVolumeID = vol->getMCid();
  119. fTrackID = gMC->GetStack()->GetCurrentTrackNumber();
  120. }
  121. // Sum energy loss for all steps in the active volume
  122. fELoss += gMC->Edep();
  123. // Create MpdStrawECTPoint at EXIT of active volume;
  124. if ( (gMC->IsTrackExiting() || gMC->IsTrackStop() || gMC->IsTrackDisappeared()) && fELoss > 0 ) {
  125. fTrackID = gMC->GetStack()->GetCurrentTrackNumber();
  126. TLorentzVector PosOut;
  127. gMC->TrackPosition(PosOut);
  128. fPosOut.SetXYZ(PosOut.X(), PosOut.Y(), PosOut.Z());
  129. //-------- define id for volumes from vol. path -------------------------------
  130. TString vol_path = gMC->CurrentVolPath();
  131. Int_t module_num = -1;
  132. Int_t submodule_num = -1;
  133. Int_t layer_num = -1;
  134. TString layer_type = "";
  135. Int_t straw_num = -1;
  136. TString key_module = "module_";
  137. TString key_submodule = "sub_module_";
  138. TString key_layer_radial = "layer_radial_";
  139. TString key_layer_stereoR = "layer_stereoR_";
  140. TString key_layer_stereoL = "layer_stereoL_";
  141. TString key_straw_tube = "straw_tube_";
  142. TString key_delim = "/";
  143. Int_t current_pos = 0;
  144. Int_t pos_in = 0;
  145. Int_t pos_out = 0;
  146. //module number
  147. pos_in = vol_path.Index(key_module, current_pos);
  148. pos_in += key_module.Length();
  149. pos_out = vol_path.Index(key_delim, pos_in);
  150. current_pos = pos_out;
  151. module_num = ((TString)(vol_path(pos_in, pos_out-pos_in))).Atoi();
  152. //submodule number
  153. pos_in = vol_path.Index(key_submodule, current_pos);
  154. pos_in += key_submodule.Length();
  155. pos_out = vol_path.Index(key_delim, pos_in);
  156. current_pos = pos_out;
  157. submodule_num = ((TString)(vol_path(pos_in, pos_out-pos_in))).Atoi();
  158. //layer number
  159. pos_in = vol_path.Index(key_layer_radial, current_pos);
  160. if(pos_in != -1) {
  161. pos_in += key_layer_radial.Length();
  162. layer_type = "radial";
  163. }
  164. if(pos_in == -1) {
  165. pos_in = vol_path.Index(key_layer_stereoR, current_pos);
  166. if(pos_in != -1) {
  167. pos_in += key_layer_stereoR.Length();
  168. layer_type = "stereoR";
  169. }
  170. }
  171. if(pos_in == -1) {
  172. pos_in = vol_path.Index(key_layer_stereoL, current_pos);
  173. if(pos_in != -1) {
  174. pos_in += key_layer_stereoL.Length();
  175. layer_type = "stereoL";
  176. }
  177. }
  178. if(pos_in != -1) {
  179. pos_out = vol_path.Index(key_delim, pos_in);
  180. current_pos = pos_out;
  181. layer_num = ((TString)(vol_path(pos_in, pos_out-pos_in))).Atoi();
  182. }
  183. //straw number
  184. pos_in = vol_path.Index(key_straw_tube, current_pos);
  185. pos_in += key_straw_tube.Length();
  186. pos_out = vol_path.Index(key_delim, pos_in);
  187. current_pos = pos_out;
  188. straw_num = ((TString)(vol_path(pos_in, pos_out-pos_in))).Atoi();
  189. //cout << "vol_path = " << vol_path << "\n";
  190. //cout << "module_num = " << module_num << "\n";
  191. //cout << "submodule_num = " << submodule_num << "\n";
  192. //cout << "layer_type = " << layer_type << "\n";
  193. //cout << "layer_num = " << layer_num << "\n";
  194. //cout << "straw_num = " << straw_num << "\n";
  195. //------------------------------------------------------------------------------
  196. MpdStrawECTPoint *p =
  197. AddHit(fTrackID, fVolumeID, fPosIn, TVector3(fMom.Px(), fMom.Py(), fMom.Pz()),
  198. fTime, (fLength+gMC->TrackLength())/2, fELoss);
  199. p->SetModule(module_num);
  200. p->SetSubmodule(submodule_num);
  201. p->SetLayer(layer_num);
  202. p->SetLayerType(layer_type);
  203. p->SetStraw(straw_num);
  204. ((MpdStack*)gMC->GetStack())->AddPoint(kECT);
  205. }
  206. return kTRUE;
  207. }
  208. //------------------------------------------------------------------------------------------------------------------------
  209. void MpdStrawECT::EndOfEvent()
  210. {
  211. if(fVerboseLevel) Print();
  212. fStrawECTCollection->Delete();
  213. fPosIndex = 0;
  214. }
  215. //------------------------------------------------------------------------------------------------------------------------
  216. void MpdStrawECT::Register(){ FairRootManager::Instance()->Register("StrawECTPoint", "StrawECT", fStrawECTCollection, kTRUE); }
  217. //------------------------------------------------------------------------------------------------------------------------
  218. TClonesArray* MpdStrawECT::GetCollection(Int_t iColl) const {
  219. if(iColl == 0) return fStrawECTCollection;
  220. return NULL;
  221. }
  222. //------------------------------------------------------------------------------------------------------------------------
  223. void MpdStrawECT::Print() const
  224. {
  225. Int_t nHits = fStrawECTCollection->GetEntriesFast();
  226. cout << "-I- MpdStrawECT: " << nHits << " points registered in this event." << endl;
  227. if(fVerboseLevel > 1)
  228. for(Int_t i=0; i<nHits; i++) (*fStrawECTCollection)[i]->Print();
  229. }
  230. //------------------------------------------------------------------------------------------------------------------------
  231. void MpdStrawECT::Reset(){ fStrawECTCollection->Delete(); ResetParameters(); }
  232. //------------------------------------------------------------------------------------------------------------------------
  233. void MpdStrawECT::CopyClones(TClonesArray* cl1, TClonesArray* cl2, Int_t offset)
  234. {
  235. Int_t nEntries = cl1->GetEntriesFast();
  236. cout << "-I- MpdStrawECT: " << nEntries << " entries to add." << endl;
  237. TClonesArray& clref = *cl2;
  238. MpdStrawECTPoint* oldpoint = NULL;
  239. for(Int_t i=0; i<nEntries; i++) {
  240. oldpoint = (MpdStrawECTPoint*) cl1->At(i);
  241. Int_t index = oldpoint->GetTrackID() + offset;
  242. oldpoint->SetTrackID(index);
  243. new (clref[fPosIndex]) MpdStrawECTPoint(*oldpoint);
  244. fPosIndex++;
  245. }
  246. cout << "-I- MpdStrawECT: " << cl2->GetEntriesFast() << " merged entries." << endl;
  247. }
  248. //------------------------------------------------------------------------------------------------------------------------
  249. void MpdStrawECT::ConstructGeometry() {
  250. TString fileName = GetGeometryFileName();
  251. if ( fileName.EndsWith(".root") ) {
  252. gLogger->Info(MESSAGE_ORIGIN,
  253. "Constructing straw ECT geometry from ROOT file %s",
  254. fileName.Data());
  255. ConstructRootGeometry();
  256. }
  257. else if ( fileName.EndsWith(".geo") ) {
  258. gLogger->Info(MESSAGE_ORIGIN,
  259. "Constructing straw ECT geometry from ASCII file %s",
  260. fileName.Data());
  261. ConstructAsciiGeometry();
  262. }
  263. else if ( fileName.EndsWith(".gdml") ) {
  264. gLogger->Info(MESSAGE_ORIGIN,
  265. "Constructing straw ECT geometry from GDML file %s",
  266. fileName.Data());
  267. //ConstructGDMLGeometry();
  268. }
  269. else {
  270. gLogger->Fatal(MESSAGE_ORIGIN,
  271. "Geometry format of straw ECT file %s not supported.",
  272. fileName.Data());
  273. }
  274. }
  275. // ----- ConstructAsciiGeometry -------------------------------------------
  276. void MpdStrawECT::ConstructAsciiGeometry() {
  277. Int_t count=0;
  278. Int_t count_tot=0;
  279. FairGeoLoader* geoLoad = FairGeoLoader::Instance();
  280. FairGeoInterface* geoFace = geoLoad->getGeoInterface();
  281. MpdStrawECTGeo* sttGeo = new MpdStrawECTGeo();
  282. sttGeo->setGeomFile(GetGeometryFileName());
  283. geoFace->addGeoModule(sttGeo);
  284. Bool_t rc = geoFace->readSet(sttGeo);
  285. if(rc) sttGeo->create(geoLoad->getGeoBuilder());
  286. TList* volList = sttGeo->getListOfVolumes();
  287. // store geo parameter
  288. FairRun *fRun = FairRun::Instance();
  289. FairRuntimeDb *rtdb = FairRun::Instance()->GetRuntimeDb();
  290. MpdStrawECTGeoPar* par =(MpdStrawECTGeoPar*)(rtdb->getContainer("MpdStrawECTGeoPar"));
  291. TObjArray *fSensNodes = par->GetGeoSensitiveNodes();
  292. TObjArray *fPassNodes = par->GetGeoPassiveNodes();
  293. FairGeoNode *node = NULL;
  294. FairGeoVolume *aVol = NULL;
  295. TListIter iter(volList);
  296. while((node = (FairGeoNode*)iter.Next())) {
  297. aVol = dynamic_cast<FairGeoVolume*> (node);
  298. if(node->isSensitive()){ fSensNodes->AddLast(aVol); count++; }
  299. else fPassNodes->AddLast(aVol);
  300. count_tot++;
  301. }
  302. par->setChanged();
  303. par->setInputVersion(fRun->GetRunId(), 1);
  304. ProcessNodes(volList);
  305. }
  306. // ----------------------------------------------------------------------------
  307. // ----- ConstructGDMLGeometry -------------------------------------------
  308. void MpdStrawECT::ConstructGDMLGeometry()
  309. {
  310. TFile *old = gFile;
  311. TGDMLParse parser;
  312. TGeoVolume* gdmlTop;
  313. // Before importing GDML
  314. Int_t maxInd = gGeoManager->GetListOfMedia()->GetEntries() - 1;
  315. gdmlTop = parser.GDMLReadFile(GetGeometryFileName());
  316. // Cheating - reassigning media indices after GDML import (need to fix this in TGDMLParse class!!!)
  317. // for (Int_t i=0; i<gGeoManager->GetListOfMedia()->GetEntries(); i++)
  318. // gGeoManager->GetListOfMedia()->At(i)->Dump();
  319. // After importing GDML
  320. Int_t j = gGeoManager->GetListOfMedia()->GetEntries() - 1;
  321. Int_t curId;
  322. TGeoMedium* m;
  323. do {
  324. m = (TGeoMedium*)gGeoManager->GetListOfMedia()->At(j);
  325. curId = m->GetId();
  326. m->SetId(curId+maxInd);
  327. j--;
  328. } while (curId > 1);
  329. // LOG(DEBUG) << "====================================================================" << FairLogger::endl;
  330. // for (Int_t i=0; i<gGeoManager->GetListOfMedia()->GetEntries(); i++)
  331. // gGeoManager->GetListOfMedia()->At(i)->Dump();
  332. Int_t newMaxInd = gGeoManager->GetListOfMedia()->GetEntries() - 1;
  333. gGeoManager->GetTopVolume()->AddNode(gdmlTop, 1, 0);
  334. ExpandNodeForGdml(gGeoManager->GetTopVolume()->GetNode(gGeoManager->GetTopVolume()->GetNdaughters()-1));
  335. for (Int_t k = maxInd+1; k < newMaxInd+1; k++) {
  336. TGeoMedium* medToDel = (TGeoMedium*)(gGeoManager->GetListOfMedia()->At(maxInd+1));
  337. LOG(DEBUG) << " removing media " << medToDel->GetName() << " with id " << medToDel->GetId() << " (k=" << k << ")" << FairLogger::endl;
  338. gGeoManager->GetListOfMedia()->Remove(medToDel);
  339. }
  340. gGeoManager->SetAllIndex();
  341. gFile = old;
  342. }
  343. void MpdStrawECT::ExpandNodeForGdml(TGeoNode* node)
  344. {
  345. LOG(DEBUG) << "----------------------------------------- ExpandNodeForGdml for node " << node->GetName() << FairLogger::endl;
  346. TGeoVolume* curVol = node->GetVolume();
  347. LOG(DEBUG) << " volume: " << curVol->GetName() << FairLogger::endl;
  348. if (curVol->IsAssembly()) {
  349. LOG(DEBUG) << " skipping volume-assembly" << FairLogger::endl;
  350. }
  351. else {
  352. TGeoMedium* curMed = curVol->GetMedium();
  353. TGeoMaterial* curMat = curVol->GetMaterial();
  354. TGeoMedium* curMedInGeoManager = gGeoManager->GetMedium(curMed->GetName());
  355. TGeoMaterial* curMatOfMedInGeoManager = curMedInGeoManager->GetMaterial();
  356. TGeoMaterial* curMatInGeoManager = gGeoManager->GetMaterial(curMat->GetName());
  357. // Current medium and material assigned to the volume from GDML
  358. LOG(DEBUG2) << " curMed\t\t\t\t" << curMed << "\t" << curMed->GetName() << "\t" << curMed->GetId() << FairLogger::endl;
  359. LOG(DEBUG2) << " curMat\t\t\t\t" << curMat << "\t" << curMat->GetName() << "\t" << curMat->GetIndex() << FairLogger::endl;
  360. // Medium and material found in the gGeoManager - either the pre-loaded one or one from GDML
  361. LOG(DEBUG2) << " curMedInGeoManager\t\t" << curMedInGeoManager
  362. << "\t" << curMedInGeoManager->GetName() << "\t" << curMedInGeoManager->GetId() << FairLogger::endl;
  363. LOG(DEBUG2) << " curMatOfMedInGeoManager\t\t" << curMatOfMedInGeoManager
  364. << "\t" << curMatOfMedInGeoManager->GetName() << "\t" << curMatOfMedInGeoManager->GetIndex() << FairLogger::endl;
  365. LOG(DEBUG2) << " curMatInGeoManager\t\t" << curMatInGeoManager
  366. << "\t" << curMatInGeoManager->GetName() << "\t" << curMatInGeoManager->GetIndex() << FairLogger::endl;
  367. TString matName = curMat->GetName();
  368. TString medName = curMed->GetName();
  369. if (curMed->GetId() != curMedInGeoManager->GetId()) {
  370. if (fFixedMedia.find(medName) == fFixedMedia.end()) {
  371. LOG(DEBUG) << " Medium needs to be fixed" << FairLogger::endl;
  372. fFixedMedia[medName] = curMedInGeoManager;
  373. Int_t ind = curMat->GetIndex();
  374. gGeoManager->RemoveMaterial(ind);
  375. LOG(DEBUG) << " removing material " << curMat->GetName()
  376. << " with index " << ind << FairLogger::endl;
  377. for (Int_t i=ind; i<gGeoManager->GetListOfMaterials()->GetEntries(); i++) {
  378. TGeoMaterial* m = (TGeoMaterial*)gGeoManager->GetListOfMaterials()->At(i);
  379. m->SetIndex(m->GetIndex()-1);
  380. }
  381. LOG(DEBUG) << " Medium fixed" << FairLogger::endl;
  382. }
  383. else {
  384. LOG(DEBUG) << " Already fixed medium found in the list " << FairLogger::endl;
  385. }
  386. }
  387. else {
  388. if (fFixedMedia.find(medName) == fFixedMedia.end()) {
  389. LOG(DEBUG) << " There is no correct medium in the memory yet" << FairLogger::endl;
  390. FairGeoLoader* geoLoad = FairGeoLoader::Instance();
  391. FairGeoInterface* geoFace = geoLoad->getGeoInterface();
  392. FairGeoMedia* geoMediaBase = geoFace->getMedia();
  393. FairGeoBuilder* geobuild = geoLoad->getGeoBuilder();
  394. FairGeoMedium* curMedInGeo = geoMediaBase->getMedium(medName);
  395. if (curMedInGeo == 0) {
  396. LOG(FATAL) << " Media not found in Geo file: " << medName << FairLogger::endl;
  397. //! This should not happen.
  398. //! This means that somebody uses material in GDML that is not in the media.geo file.
  399. //! Most probably this is the sign to the user to check materials' names in the CATIA model.
  400. }
  401. else {
  402. LOG(DEBUG) << " Found media in Geo file" << medName << FairLogger::endl;
  403. Int_t nmed = geobuild->createMedium(curMedInGeo);
  404. fFixedMedia[medName] = (TGeoMedium*)gGeoManager->GetListOfMedia()->Last();
  405. gGeoManager->RemoveMaterial(curMatOfMedInGeoManager->GetIndex());
  406. LOG(DEBUG) << " removing material " << curMatOfMedInGeoManager->GetName()
  407. << " with index " << curMatOfMedInGeoManager->GetIndex() << FairLogger::endl;
  408. for (Int_t i=curMatOfMedInGeoManager->GetIndex(); i<gGeoManager->GetListOfMaterials()->GetEntries(); i++) {
  409. TGeoMaterial* m = (TGeoMaterial*)gGeoManager->GetListOfMaterials()->At(i);
  410. m->SetIndex(m->GetIndex()-1);
  411. }
  412. }
  413. if (curMedInGeo->getSensitivityFlag()) {
  414. LOG(DEBUG) << " Adding sensitive " << curVol->GetName() << FairLogger::endl;
  415. AddSensitiveVolume(curVol);
  416. }
  417. }
  418. else {
  419. LOG(DEBUG) << " Already fixed medium found in the list" << FairLogger::endl;
  420. LOG(DEBUG) << "!!! Sensitivity: " << fFixedMedia[medName]->GetParam(0) << FairLogger::endl;
  421. if (fFixedMedia[medName]->GetParam(0) == 1) {
  422. LOG(DEBUG) << " Adding sensitive " << curVol->GetName() << FairLogger::endl;
  423. AddSensitiveVolume(curVol);
  424. }
  425. }
  426. }
  427. curVol->SetMedium(fFixedMedia[medName]);
  428. gGeoManager->SetAllIndex();
  429. // gGeoManager->GetListOfMaterials()->Print();
  430. // gGeoManager->GetListOfMedia()->Print();
  431. }
  432. //! Recursevly go down the tree of nodes
  433. if (curVol->GetNdaughters() != 0) {
  434. TObjArray* NodeChildList = curVol->GetNodes();
  435. TGeoNode* curNodeChild;
  436. for (Int_t j=0; j<NodeChildList->GetEntriesFast(); j++) {
  437. curNodeChild = (TGeoNode*)NodeChildList->At(j);
  438. ExpandNodeForGdml(curNodeChild);
  439. }
  440. }
  441. }
  442. //-----------------------------------------------------------------------------
  443. //Check if Sensitive-----------------------------------------------------------
  444. Bool_t MpdStrawECT::CheckIfSensitive(std::string name) {
  445. TString tsname = name;
  446. if (tsname.Contains("straw_gas")) {
  447. return kTRUE;
  448. }
  449. return kFALSE;
  450. }
  451. //-----------------------------------------------------------------------------
  452. MpdStrawECTPoint* MpdStrawECT::AddHit(Int_t trackID, Int_t detID, TVector3 pos,
  453. TVector3 mom, Double_t time, Double_t length, Double_t eLoss) {
  454. TClonesArray& clref = *fStrawECTCollection;
  455. Int_t size = clref.GetEntriesFast();
  456. //std::cout << "ELoss: " << eLoss << "\n";
  457. return new(clref[size]) MpdStrawECTPoint(trackID, detID, pos, mom, time, length, eLoss);
  458. }
  459. //------------------------------------------------------------------------------------------------------------------------
  460. ClassImp(MpdStrawECT)