#include <Geometry.h>
Public Member Functions | |
const PlaneGeo & | Plane (unsigned int i) const |
const std::set< unsigned int > & | GetPlanesByView (View_t v=kXorY) const |
const unsigned int | NextPlaneInView (unsigned int p, int d=+1) const |
const unsigned int | NextPlaneOtherView (unsigned int p, int d=+1) const |
const CellGeo & | IdToCell (const CellUniqueId &id, int *iplane=0, int *icell=0) const |
double | DetHalfWidth () const |
double | DetHalfHeight () const |
double | DetLength () const |
double | SurfaceY () const |
void | WorldBox (double *xlo, double *xhi, double *ylo, double *yhi, double *zlo, double *zhi) const |
int | CurrentCell (int *ip, int *ic) const |
const CellUniqueId | CurrentCellId () const |
const TGeoMaterial * | Material (double x, double y, double z) const |
double | TotalMass () const |
double | MassBetweenPoints (double *p1, double *p2) const |
TGeoManager * | ROOTGeoManager () const |
Static Public Member Functions | |
Geometry & | Instance (int runnum, short detId) |
Private Types | |
typedef std::pair< unsigned short, unsigned short > | PlaneCellPair |
typedef std::map< CellUniqueId, PlaneCellPair > | UIDMap |
Private Member Functions | |
Geometry (const char *gdmlfile) | |
~Geometry () | |
void | FindPlanes (std::vector< const TGeoNode * > &n, unsigned int depth) |
void | MakePlane (std::vector< const TGeoNode * > &n, unsigned int depth) |
void | SetDrawOptions () |
void | BuildMaps () |
void | SetDetectorSize () |
Private Attributes | |
int | fRunNumber |
Run number of configuration. | |
int | fDetId |
Id number of detector geometry. | |
std::vector< PlaneGeo * > | fPlanes |
The detector planes. | |
UIDMap | fIdMap |
Unique ID -> Plane,Cell. | |
std::set< unsigned int > | fAllPlanes |
List of all planes. | |
std::set< unsigned int > | fXplanes |
List of X measuring planes. | |
std::set< unsigned int > | fYplanes |
List of Y measuring planes. | |
double | fDetLength |
Detector 1/2 length (cm). | |
double | fDetHalfHeight |
Detector 1/2 height (cm). | |
double | fDetHalfWidth |
Detector 1/2 width (cm). | |
Static Private Attributes | |
Geometry * | fInstance = 0 |
Instance of geometry. |
Definition at line 30 of file Geometry.h.
|
Definition at line 73 of file Geometry.h. Referenced by BuildMaps(). |
|
Definition at line 74 of file Geometry.h. |
|
Definition at line 91 of file Geometry.cxx. References BuildMaps(), FindPlanes(), fPlanes, SetDetectorSize(), and SetDrawOptions(). Referenced by Instance(). 00092 { 00093 std::string fname; 00094 00095 // Test if the gdml file can be openned using exactly the name we've 00096 // been passed 00097 while (1) { 00098 struct stat sb; 00099 fname = gdmlfile; 00100 if (stat(fname.c_str(), &sb)==0) break; 00101 00102 const char* srt_private = getenv("SRT_PRIVATE_CONTEXT"); 00103 if (srt_private!=0) { 00104 fname = srt_private; 00105 fname += "/Geometry/gdml/"; 00106 fname += gdmlfile; 00107 if (stat(fname.c_str(), &sb)==0) break; 00108 } 00109 00110 const char* srt_public = getenv("SRT_PUBLIC_CONTEXT"); 00111 if (srt_public!=0) { 00112 fname = srt_public; 00113 fname += "/Geometry/gdml/"; 00114 fname += gdmlfile; 00115 if (stat(fname.c_str(), &sb)==0) break; 00116 } 00117 00118 // Failed to resolve the file name 00119 throw Exception(__FILE__, __LINE__, Exception::FILE_NOT_FOUND, 00120 "GDML file not found."); 00121 } 00122 00123 // Open the GDML file 00124 TGeoManager::Import(fname.c_str()); 00125 this->SetDrawOptions(); 00126 00127 std::vector<const TGeoNode*> n(16); 00128 n[0] = gGeoManager->GetTopNode(); 00129 this->FindPlanes(n, 0); 00130 sort(fPlanes.begin(), fPlanes.end(), plane_sort); 00131 00132 this->SetDetectorSize(); 00133 this->BuildMaps(); 00134 }
|
|
Definition at line 144 of file Geometry.cxx. References fPlanes. 00145 { 00146 for (unsigned int i=0; i<fPlanes.size(); ++i) { 00147 if (fPlanes[i]) { delete fPlanes[i]; fPlanes[i] = 0; } 00148 } 00149 }
|
|
Build maps used for quick access to the geometry Definition at line 424 of file Geometry.cxx. References geo::PlaneGeo::Cell(), fAllPlanes, fIdMap, fPlanes, fXplanes, fYplanes, geo::CellGeo::Id(), geo::PlaneGeo::Ncells(), Plane(), and PlaneCellPair. Referenced by Geometry(). 00425 { 00426 // Sets which list planes by view 00427 for (unsigned int i=0; i<fPlanes.size(); ++i) { 00428 fAllPlanes.insert(i); 00429 if (fPlanes[i]->View()==kX) fXplanes.insert(i); 00430 else if (fPlanes[i]->View()==kY) fYplanes.insert(i); 00431 else { 00432 throw Exception(__FILE__,__LINE__,Exception::BAD_GEO_CONFIG, 00433 "Bad plane view!"); 00434 } 00435 } 00436 00437 // Unique Id -> cell map 00438 for (unsigned int i=0; i<fPlanes.size(); ++i) { 00439 const PlaneGeo* p = fPlanes[i]; 00440 for (int j=0; j<p->Ncells(); ++j) { 00441 // Assert that the cell ids are in fact unique 00442 if (fIdMap.find(this->Plane(i).Cell(j).Id())!=fIdMap.end()) { 00443 throw Exception(__FILE__,__LINE__,Exception::BAD_CELL_UNIQUE_ID, 00444 "Duplicate Cell UIDs"); 00445 } 00446 PlaneCellPair p(i,j); 00447 fIdMap[this->Plane(i).Cell(j).Id()] = p; 00448 } 00449 } 00450 }
|
|
Return the plane and cell number for the current tracking position held by the gGeoManager
Definition at line 160 of file Geometry.cxx. References IdToCell(). Referenced by testFindCell(). 00161 { 00162 try { 00163 this->IdToCell(this->CurrentCellId(), planeid, cellid); 00164 } 00165 catch (Exception e) { 00166 *planeid = -1; 00167 *cellid = -1; 00168 return -1; 00169 } 00170 return 1; 00171 }
|
|
Return the unique cell Id# (well, #'s) for the current cell Definition at line 177 of file Geometry.cxx. References geo::CellUniqueId, and geo::PathToUniqueId(). 00178 { 00179 return PathToUniqueId(gGeoManager->GetPath()); 00180 }
|
|
Definition at line 199 of file Geometry.cxx. Referenced by evd::DetectorView::DetectorView(), evd::Display3D::Draw(), evgen::GENIEHelper::InitializeFluxDriver(), main(), and evd::PlaneView::PlaneView(). 00200 {
00201 return fDetHalfHeight;
00202 }
|
|
Definition at line 192 of file Geometry.cxx. Referenced by evd::DetectorView::DetectorView(), evd::Display3D::Draw(), main(), and evd::PlaneView::PlaneView(). 00193 {
00194 return fDetHalfWidth;
00195 }
|
|
Definition at line 206 of file Geometry.cxx. Referenced by mcchk::CosmicAna::Ana(), evd::DetectorView::DetectorView(), evd::Display3D::Draw(), evgen::GENIEHelper::InitializeFluxDriver(), main(), evgen::GENIEHelper::PackMCTruth(), evd::PlaneView::PlaneView(), evgen::CRYGen::Sample(), and evgen::GENIEHelper::Update(). 00207 {
00208 return fDetLength;
00209 }
|
|
Recursively search for planes in the geometry
Definition at line 281 of file Geometry.cxx. References MakePlane(). Referenced by Geometry(). 00283 { 00284 const char* nm = n[depth]->GetName(); 00285 if (strncmp(nm,"vPlane", 6)==0) { 00286 this->MakePlane(n, depth); 00287 return; 00288 } 00289 00290 // Explore the next layer down 00291 unsigned int deeper = depth+1; 00292 if (deeper>=n.size()) { 00293 throw Exception(__FILE__,__LINE__,Exception::BAD_NODE, 00294 "Exceeded maximum depth."); 00295 } 00296 const TGeoVolume* v = n[depth]->GetVolume(); 00297 int nd = v->GetNdaughters(); 00298 for (int i=0; i<nd; ++i) { 00299 n[deeper] = v->GetNode(i); 00300 this->FindPlanes(n, deeper); 00301 } 00302 }
|
|
Return list of planes which measure the requested projection
Definition at line 361 of file Geometry.cxx. Referenced by mcchk::CosmicAna::Ana(), vali::Validator::BookPT(), subshower::RecoSubShower::FindMaxWindow(), and testFindPlanes(). 00362 { 00363 // Return the user's choice, fall back on all planes 00364 if (v==kX) return fXplanes; 00365 else if (v==kY) return fYplanes; 00366 return fAllPlanes; 00367 }
|
|
Given a unique cell identifier, look up the plane and cell numbers.
Definition at line 341 of file Geometry.cxx. References geo::PlaneGeo::Cell(), fIdMap, and Plane(). Referenced by CurrentCell(), novamc::MCApplication::Stepping(), and testUniqueId(). 00344 { 00345 UIDMap::const_iterator itr = fIdMap.find(id); 00346 if (itr == fIdMap.end()) { 00347 throw Exception(__FILE__,__LINE__,Exception::BAD_CELL_UNIQUE_ID); 00348 } 00349 *iplane = itr->second.first; 00350 *icell = itr->second.second; 00351 return this->Plane(*iplane).Cell(*icell); 00352 }
|
|
Build a detector geometry
Definition at line 49 of file Geometry.cxx. References fDetId, fInstance, fRunNumber, and Geometry(). Referenced by fillntup::FillNtup::Ana(), mcchk::CosmicAna::Ana(), evd::DetectorView::DetectorView(), evd::Display3D::Draw(), evd::DetectorView::Draw(), evd::PlaneView::Draw3DShowers(), evd::PlaneView::DrawMCHits(), evd::Display3D::DrawMCHits(), evd::PlaneView::DrawPlaneClusters(), evd::Display3D::DrawRawData(), evd::PlaneView::DrawRawDigit(), evd::PlaneView::DrawShowers(), evd::PlaneView::DrawTrackCandidates(), subshower::RecoSubShower::FormStripList(), evgen::GENIEHelper::InitializeFluxDriver(), evgen::GENIEHelper::InitializeGeometry(), main(), novamc::MCApplication::MCApplication(), evgen::GENIEHelper::PackMCTruth(), evd::PlaneView::PlaneView(), calhit::CalHit::Raw(), vali::Validator::Reco(), rpr::TrackReco::Reco(), spider::SpiderWeb::Reco(), photrans::SimpleTransport::Reco(), photrans::SimpleEnhancedTransport::Reco(), subshower::RecoSubShower::Reco(), rpr::FindTrackSeg::Reco(), ctrk::CosmicTrack::Reco(), cluster::Clusterer::Reco(), and evgen::GENIEHelper::Update(). 00050 { 00051 int oldrunnum = -1; 00052 int olddetid = -1; 00053 // If our current geometry doesn't match request, get rid of it 00054 if (fInstance!=0) { 00055 if (fInstance->fRunNumber != runnum || fInstance->fDetId != detId) { 00056 oldrunnum = fInstance->fRunNumber; 00057 olddetid = fInstance->fDetId; 00058 delete fInstance; 00059 fInstance = 0; 00060 gGeoManager->Clear(); 00061 } 00062 } 00063 00064 // Build a new geometry if we need to 00065 if (fInstance==0) { 00066 static char* fname[3] = {"neardet.gdml","fardet.gdml","ipnd.gdml"}; 00067 if (detId!=rawdata::kNear && 00068 detId!=rawdata::kFar && 00069 detId!=rawdata::kIPND) { 00070 char msg[256]; 00071 sprintf(msg,"Unknown detector id number: %d",(int)detId); 00072 throw Exception(__FILE__,__LINE__,Exception::BAD_GEO_CONFIG,msg); 00073 } 00074 if (runnum<0) { 00075 char msg[256]; 00076 std::sprintf(msg,"Bad run number: %d",runnum); 00077 throw Exception(__FILE__,__LINE__,Exception::BAD_GEO_CONFIG,msg); 00078 } 00079 fInstance = new Geometry(fname[detId]); 00080 std::cout << "Geometry updated from " 00081 << olddetid << ":" << oldrunnum << " to " 00082 << detId << ":" << runnum << std::endl; 00083 fInstance->fDetId = detId; 00084 fInstance->fRunNumber = runnum; 00085 } 00086 return *fInstance; 00087 }
|
|
Definition at line 306 of file Geometry.cxx. References fPlanes. Referenced by FindPlanes().
|
|
Return the column density between 2 points
The purpose of this method is to determine the column density between the two points given. Do that by starting at p1 and stepping until you get to the node of p2. calculate the distance between the point just inside that node and p2 to get the last bit of column density first initialize a track - get the direction cosines might be helpful to have a point to a TGeoNode check that the points are not in the same volume already. if they are in different volumes, keep stepping until you are in the same volume as the second point the act of stepping puts you in the next node and returns that node now you are in the same volume as the last point, but not at that point. get the distance between the current point and the last one Definition at line 487 of file Geometry.cxx. 00488 { 00489 00495 double columnD = 0.; 00496 00498 double length = TMath::Sqrt(TMath::Power(p2[0]-p1[0], 2.) 00499 + TMath::Power(p2[1]-p1[1], 2.) 00500 + TMath::Power(p2[2]-p1[2], 2.)); 00501 double dxyz[3] = {(p2[0]-p1[0])/length, (p2[1]-p1[1])/length, (p2[2]-p1[2])/length}; 00502 00503 gGeoManager->InitTrack(p1,dxyz); 00504 00506 TGeoNode *node = gGeoManager->GetCurrentNode(); 00507 00511 while(!gGeoManager->IsSameLocation(p2[0], p2[1], p2[2])){ 00512 gGeoManager->FindNextBoundary(); 00513 columnD += gGeoManager->GetStep()*node->GetMedium()->GetMaterial()->GetDensity(); 00514 00516 node = gGeoManager->Step(); 00517 }//end loop to get to volume of second point 00518 00521 const double *current = gGeoManager->GetCurrentPoint(); 00522 length = TMath::Sqrt(TMath::Power(p2[0]-current[0], 2.) 00523 + TMath::Power(p2[1]-current[1], 2.) 00524 + TMath::Power(p2[2]-current[2], 2.)); 00525 columnD += length*node->GetMedium()->GetMaterial()->GetDensity(); 00526 00527 return columnD; 00528 }
|
|
Definition at line 184 of file Geometry.cxx. Referenced by proto_evgen::TargetVolume::RandomVertex(). 00185 { 00186 const TGeoNode* node = gGeoManager->FindNode(x,y,z); 00187 return node->GetMedium()->GetMaterial(); 00188 }
|
|
Return the index of the next plane in a particular view
Definition at line 380 of file Geometry.cxx. References fPlanes. Referenced by subshower::RecoSubShower::LinkedBlobs(), calhit::CalHit::Neighborinos(), testFindPlanes(), and subshower::RecoSubShower::TransCluster(). 00381 { 00382 assert(p1<fPlanes.size()); 00383 00384 // March along from p1 in the direction d until we find a match 00385 int s = (d>0 ? 1 : -1); 00386 if (s<0 && p1==0) return kPLANE_NOT_FOUND; 00387 if (s>0 && p1==fPlanes.size()-1) return kPLANE_NOT_FOUND; 00388 for (unsigned int p2=p1+d; 1; p2+=s) { 00389 if (fPlanes[p1]->View() == fPlanes[p2]->View()) return p2; 00390 if (p2==fPlanes.size()-1) return kPLANE_NOT_FOUND; 00391 if (p2==0) return kPLANE_NOT_FOUND; 00392 } 00393 }
|
|
Return the index of the next plane in a particular view
Definition at line 405 of file Geometry.cxx. References fPlanes. Referenced by calhit::CalHit::Neighborinos(), and testFindPlanes(). 00406 { 00407 assert(p1<fPlanes.size()); 00408 00409 // March along from p1 in the direction d until we find a match 00410 int s = (d>0 ? 1 : -1); 00411 if (s<0 && p1==0) return kPLANE_NOT_FOUND; 00412 if (s>0 && p1==fPlanes.size()-1) return kPLANE_NOT_FOUND; 00413 for (unsigned int p2=p1+d; 1; p2+=s) { 00414 if (fPlanes[p1]->View() != fPlanes[p2]->View()) return p2; 00415 if (p2==fPlanes.size()) return kPLANE_NOT_FOUND; 00416 if (p2==0) return kPLANE_NOT_FOUND; 00417 } 00418 }
|
|
|
Definition at line 137 of file Geometry.cxx. Referenced by evgen::GENIEHelper::InitializeGeometry(). 00138 {
00139 return gGeoManager;
00140 }
|
|
Calcualte the detector size based on the geomtry. Requires a little bit of work, so do this once when the geometry is updated Definition at line 536 of file Geometry.cxx. References geo::PlaneGeo::Cell(), fDetHalfHeight, fDetHalfWidth, fDetLength, fPlanes, geo::CellGeo::GetCenter(), and Plane(). Referenced by Geometry(). 00537 { 00538 double dummy; 00539 00540 // Compute the height and width based on the sizes of the vertical 00541 // and horizontal planes. Remember, in the local plane frame z goes 00542 // along the length of the modules. Hence the use of axis=3 to get 00543 // the detector height (from the vertical planes) and width (from 00544 // the horizontal planes). 00545 const TGeoVolume* vvp = gGeoManager->GetVolume("vPlaneV"); 00546 const TGeoVolume* hvp = gGeoManager->GetVolume("vPlaneH"); 00547 00548 // Handle the case of the IPND where the volumes are names slightly 00549 // different accounting for the different widths of the modules. 00550 if (vvp==0) vvp = gGeoManager->GetVolume("vPlaneVIPND"); 00551 if (hvp==0) hvp = gGeoManager->GetVolume("vPlaneHIPND"); 00552 if (vvp==0 && hvp==0) { 00553 vvp = gGeoManager->GetVolume("vPlaneVIPND"); 00554 hvp = gGeoManager->GetVolume("vPlaneHIPND"); 00555 } 00556 if (vvp==0 || hvp==0) { 00557 throw Exception(__FILE__,__LINE__, 00558 Exception::BAD_GEO_CONFIG, 00559 "Unable to find shapes to set detector size"); 00560 } 00561 00562 const TGeoShape* vp = vvp->GetShape(); 00563 const TGeoShape* hp = hvp->GetShape(); 00564 vp->GetAxisRange(3,dummy,fDetHalfHeight); 00565 hp->GetAxisRange(3,dummy,fDetHalfWidth); 00566 00567 // Remember -- in the local plane coordinate frame, the "depth" is 00568 // along y, hence axis=2. In the world coordinate this is along z. 00569 double vplanez1, vplanez2; 00570 vp->GetAxisRange(2,vplanez1,vplanez2); 00571 00572 // Compute the detector length based on the center positions of 00573 // cells in the first and last plane. Adjust to account for the 00574 // plane widths. This ajustment works if the vertical and horizontal 00575 // planes have the same depth, or if the detector begins and ends 00576 // with vertical planes. As of the TDR both of these conditions was 00577 // true for all the NOvA detectors. 00578 double xyz1[3] = {0,0,0}; 00579 double xyz2[3] = {0,0,0}; 00580 this->Plane(0). Cell(0).GetCenter(xyz1); 00581 this->Plane(fPlanes.size()-1).Cell(0).GetCenter(xyz2); 00582 fDetLength = (xyz2[2]-xyz1[2])+(vplanez2-vplanez1); 00583 }
|
|
Definition at line 270 of file Geometry.cxx. Referenced by Geometry(). 00270 { }
|
|
A typical y-position value at the surface (where earth meets air) for this detector site
Definition at line 224 of file Geometry.cxx. Referenced by evgen::CRYGen::Sample(). 00225 { 00226 switch (fDetId) { 00227 case 0: return 130.0E2; // near det: surface ~130m up from det. center 00228 case 1: return 9.94E2; // far det: surface ~10m up from det. center 00229 case 2: return -1.966E2; // IPND: surface ~2 meter below det, center 00230 default: throw Exception(__FILE__,__LINE__,Exception::BAD_GEO_CONFIG); 00231 } 00232 return 0.0; 00233 }
|
|
Return the total mass of the detector the detector resides in the vDetEnclosure volume so start there loop over the daughters in the volume and let ROOT calculate the mass in kg for you for the daughter nodes that have vBlock in their names this only goes one level down, so no worries of double counting non-superblock blocks Definition at line 457 of file Geometry.cxx. Referenced by main(), and evgen::GENIEHelper::Update(). 00458 { 00459 double mass = 0.; 00460 00462 const TGeoVolume *enclosure = gGeoManager->FindVolumeFast("vDetEnclosure"); 00463 00468 for(int d = 0; d < enclosure->GetNdaughters(); ++d){ 00469 const TGeoNode* node = enclosure->GetNode(d); 00470 const char* nm = node->GetName(); 00471 if(strncmp(nm,"vSuperBlock",11) == 0){ 00472 std::string name = nm; 00473 name.erase(name.find('_')); 00474 mass += gGeoManager->FindVolumeFast(name.c_str())->Weight(); 00475 } 00476 } 00477 return mass; 00478 }
|
|
Return the ranges of x,y and z for the "world volume" that the entire geometry lives in. If any pointers are 0, then those coordinates are ignored.
Definition at line 247 of file Geometry.cxx. Referenced by evgen::CRYGen::Sample(), and testProject(). 00250 { 00251 const TGeoShape* s = gGeoManager->GetVolume("vWorld")->GetShape(); 00252 assert(s); 00253 00254 if (xlo || xhi) { 00255 double x1, x2; 00256 s->GetAxisRange(1,x1,x2); if (xlo) *xlo = x1; if (xhi) *xhi = x2; 00257 } 00258 if (ylo || yhi) { 00259 double y1, y2; 00260 s->GetAxisRange(2,y1,y2); if (ylo) *ylo = y1; if (yhi) *yhi = y2; 00261 } 00262 if (zlo || zhi) { 00263 double z1, z2; 00264 s->GetAxisRange(3,z1,z2); if (zlo) *zlo = z1; if (zhi) *zhi = z2; 00265 } 00266 }
|
|
List of all planes.
Definition at line 82 of file Geometry.h. Referenced by BuildMaps(). |
|
Detector 1/2 height (cm).
Definition at line 87 of file Geometry.h. Referenced by SetDetectorSize(). |
|
Detector 1/2 width (cm).
Definition at line 88 of file Geometry.h. Referenced by SetDetectorSize(). |
|
Id number of detector geometry.
Definition at line 78 of file Geometry.h. Referenced by Instance(). |
|
Detector 1/2 length (cm).
Definition at line 86 of file Geometry.h. Referenced by SetDetectorSize(). |
|
Unique ID -> Plane,Cell.
Definition at line 80 of file Geometry.h. Referenced by BuildMaps(), and IdToCell(). |
|
Instance of geometry.
Definition at line 28 of file Geometry.cxx. Referenced by Instance(). |
|
The detector planes.
Definition at line 79 of file Geometry.h. Referenced by BuildMaps(), Geometry(), MakePlane(), NextPlaneInView(), NextPlaneOtherView(), Plane(), SetDetectorSize(), and ~Geometry(). |
|
Run number of configuration.
Definition at line 77 of file Geometry.h. Referenced by Instance(). |
|
List of X measuring planes.
Definition at line 83 of file Geometry.h. Referenced by BuildMaps(). |
|
List of Y measuring planes.
Definition at line 84 of file Geometry.h. Referenced by BuildMaps(). |