Main Page | Namespace List | Class Hierarchy | Alphabetical List | Class List | Directories | File List | Namespace Members | Class Members | File Members | Related Pages

geo::Geometry Class Reference

The geometry of one entire detector (near, far, ipnd). More...

#include <Geometry.h>

List of all members.

Public Member Functions

const PlaneGeoPlane (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
void CellInfo (unsigned int ip, unsigned int ic, View_t *view=0, double *pos=0, double *dpos=0) const
const CellGeoIdToCell (const CellUniqueId &id, int *iplane=0, int *icell=0) const
unsigned int NPlanes () 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

GeometryInstance (int runnum, const char *fname)

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.
std::string fGDMLFile
 gdml file holding the 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

GeometryfInstance = 0
 Instance of geometry.


Detailed Description

The geometry of one entire detector (near, far, ipnd).

Definition at line 31 of file Geometry.h.


Member Typedef Documentation

typedef std::pair<unsigned short, unsigned short> geo::Geometry::PlaneCellPair [private]
 

Definition at line 77 of file Geometry.h.

Referenced by BuildMaps().

typedef std::map<CellUniqueId, PlaneCellPair> geo::Geometry::UIDMap [private]
 

Definition at line 78 of file Geometry.h.


Constructor & Destructor Documentation

Geometry::Geometry const char *  gdmlfile  )  [private]
 

Definition at line 86 of file Geometry.cxx.

References BuildMaps(), FindPlanes(), fPlanes, SetDetectorSize(), and SetDrawOptions().

Referenced by Instance().

00087 {
00088   std::string fname;
00089 
00090   // Test if the gdml file can be openned using exactly the name we've
00091   // been passed
00092   while (1) {
00093     struct stat sb;
00094     fname = gdmlfile;
00095     if (stat(fname.c_str(), &sb)==0) break;
00096 
00097     const char* srt_private = getenv("SRT_PRIVATE_CONTEXT");
00098     if (srt_private!=0) {
00099       fname = srt_private;
00100       fname += "/Geometry/gdml/";
00101       fname += gdmlfile;
00102       if (stat(fname.c_str(), &sb)==0) break;
00103     }
00104 
00105     const char* srt_public = getenv("SRT_PUBLIC_CONTEXT");
00106     if (srt_public!=0) {
00107       fname = srt_public;
00108       fname += "/Geometry/gdml/";
00109       fname += gdmlfile;
00110       if (stat(fname.c_str(), &sb)==0) break;
00111     }
00112     
00113     // Failed to resolve the file name
00114     throw Exception(__FILE__, __LINE__, Exception::FILE_NOT_FOUND,
00115                     "GDML file not found.");
00116   }
00117 
00118   // Open the GDML file
00119   TGeoManager::Import(fname.c_str());
00120   this->SetDrawOptions();
00121   
00122   std::vector<const TGeoNode*> n(16);
00123   n[0] = gGeoManager->GetTopNode();
00124   this->FindPlanes(n, 0);
00125   sort(fPlanes.begin(), fPlanes.end(), plane_sort);
00126   
00127   this->SetDetectorSize();
00128   this->BuildMaps();
00129 
00130   std::cout << "Detector geometry loaded from " << fname << std::endl;
00131 }

Geometry::~Geometry  )  [private]
 

Definition at line 142 of file Geometry.cxx.

References fPlanes.

00143 {
00144   for (unsigned int i=0; i<fPlanes.size(); ++i) {
00145     if (fPlanes[i]) { delete fPlanes[i]; fPlanes[i] = 0; }
00146   }
00147 }


Member Function Documentation

void Geometry::BuildMaps  )  [private]
 

Build maps used for quick access to the geometry

Definition at line 463 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().

00464 {
00465   // Sets which list planes by view
00466   for (unsigned int i=0; i<fPlanes.size(); ++i) {
00467     fAllPlanes.insert(i);
00468     if      (fPlanes[i]->View()==kX) fXplanes.insert(i);
00469     else if (fPlanes[i]->View()==kY) fYplanes.insert(i);
00470     else {
00471       throw Exception(__FILE__,__LINE__,Exception::BAD_GEO_CONFIG,
00472                       "Bad plane view!");
00473     }
00474   }
00475   
00476   // Unique Id -> cell map
00477   for (unsigned int i=0; i<fPlanes.size(); ++i) {
00478     const PlaneGeo* p = fPlanes[i];
00479     for (unsigned int j=0; j<p->Ncells(); ++j) {
00480       // Assert that the cell ids are in fact unique
00481       if (fIdMap.find(this->Plane(i).Cell(j).Id())!=fIdMap.end()) {
00482         throw Exception(__FILE__,__LINE__,Exception::BAD_CELL_UNIQUE_ID,
00483                         "Duplicate Cell UIDs");
00484       }
00485       PlaneCellPair p(i,j);
00486       fIdMap[this->Plane(i).Cell(j).Id()] = p;
00487     }
00488   }
00489 }

void Geometry::CellInfo unsigned int  ip,
unsigned int  ic,
View_t view = 0,
double *  pos = 0,
double *  dpos = 0
const
 

Given a plane and cell number, return which view the cell measures, its nominal center position, and nominal size about that center position.

Parameters:
ip - plane number
ic - cell number
view - on retrun, the view this cell measures in
pos - the nominal (xyz) center position of this cell (cm)
dpos - half-dimensions of the cell in xyz (cm)
Exceptions:
geo::Exception 

Definition at line 339 of file Geometry.cxx.

References geo::PlaneGeo::Cell(), geo::CellGeo::GetCenter(), geo::CellGeo::HalfD(), geo::CellGeo::HalfL(), geo::CellGeo::HalfW(), geo::PlaneGeo::Ncells(), NPlanes(), Plane(), and geo::PlaneGeo::View().

Referenced by evd::Display3D::DrawRawData(), and evd::RecoBaseDrawer::GetClusterOutlines().

00341 {
00342   if (ip<0 || ip>=this->NPlanes()) {
00343     throw Exception(__FILE__,__LINE__,Exception::BAD_PLANE_NUMBER);
00344   }
00345   const PlaneGeo& p = this->Plane(ip);
00346   if (view) *view = p.View();
00347   
00348   if (ic<0 || ic>=p.Ncells()) {
00349     throw Exception(__FILE__,__LINE__,Exception::BAD_CELL_NUMBER);
00350   }
00351   const::CellGeo& c = p.Cell(ic);
00352   if (pos) c.GetCenter(pos);
00353   
00354   if (dpos) {
00355     if (p.View() == kX) {
00356       dpos[0] = c.HalfW();
00357       dpos[1] = c.HalfL();
00358       dpos[2] = c.HalfD();
00359     }
00360     else if (p.View() == kY) {
00361       dpos[0] = c.HalfL();
00362       dpos[1] = c.HalfW();
00363       dpos[2] = c.HalfD();
00364     }
00365   }
00366 }

int Geometry::CurrentCell int *  planeid,
int *  cellid
const
 

Return the plane and cell number for the current tracking position held by the gGeoManager

Parameters:
planeid : pointer to plane number, >=0 on success, -1 on failure
cellid : pointer to cell number, >=0 on success, -1 on failure
Returns:
>0 if cell is found, <0 if not found.

Definition at line 158 of file Geometry.cxx.

References IdToCell().

Referenced by testFindCell().

00159 {
00160   try {
00161     this->IdToCell(this->CurrentCellId(), planeid, cellid);
00162   }
00163   catch (Exception e) {
00164     *planeid = -1;
00165     *cellid  = -1;
00166     return -1;
00167   }
00168   return 1;
00169 }

const CellUniqueId Geometry::CurrentCellId  )  const
 

Return the unique cell Id# (well, #'s) for the current cell

Definition at line 175 of file Geometry.cxx.

References geo::CellUniqueId, and geo::PathToUniqueId().

00176 {
00177   return PathToUniqueId(gGeoManager->GetPath());
00178 }

double Geometry::DetHalfHeight  )  const
 

Definition at line 197 of file Geometry.cxx.

Referenced by mcchk::NeutrinoAna::Ana(), evd::DetectorView::DetectorView(), evd::Display3D::Draw(), main(), and evd::TZProjPad::TZProjPad().

00198 {
00199   return fDetHalfHeight;
00200 }

double Geometry::DetHalfWidth  )  const
 

Definition at line 190 of file Geometry.cxx.

Referenced by mcchk::NeutrinoAna::Ana(), evd::DetectorView::DetectorView(), evd::Display3D::Draw(), main(), and evd::TZProjPad::TZProjPad().

00191 {
00192   return fDetHalfWidth;
00193 }

double Geometry::DetLength  )  const
 

Definition at line 204 of file Geometry.cxx.

Referenced by mcchk::NeutrinoAna::Ana(), mcchk::CosmicAna::Ana(), clusterss::CompareClusters::Ana(), clusterss::ClusterCheck::Ana(), evd::DetectorView::DetectorView(), evd::Display3D::Draw(), and evd::TZProjPad::TZProjPad().

00205 { 
00206   return fDetLength;
00207 }

void Geometry::FindPlanes std::vector< const TGeoNode * > &  n,
unsigned int  depth
[private]
 

Recursively search for planes in the geometry

Parameters:
n - Array holding the path of nodes tranversed
depth - Depth of the search to this point
Exceptions:
geo::Exception is maximum depth is exceeded.

Definition at line 278 of file Geometry.cxx.

References MakePlane().

Referenced by Geometry().

00280 {
00281   const char* nm = n[depth]->GetName();
00282   if (strncmp(nm,"vPlane", 6)==0) {
00283     this->MakePlane(n, depth);
00284     return;
00285   }
00286 
00287   // Explore the next layer down
00288   unsigned int deeper = depth+1;
00289   if (deeper>=n.size()) {
00290     throw Exception(__FILE__,__LINE__,Exception::BAD_NODE,
00291                     "Exceeded maximum depth.");
00292   }
00293   const TGeoVolume* v = n[depth]->GetVolume();
00294   int nd = v->GetNdaughters();
00295   for (int i=0; i<nd; ++i) {
00296     n[deeper] = v->GetNode(i);
00297     this->FindPlanes(n, deeper);
00298   }
00299 }

const std::set< unsigned int > & Geometry::GetPlanesByView View_t  v = kXorY  )  const
 

Return list of planes which measure the requested projection

Parameters:
v : X/Y or All (see View_t typedef in PlaneGeo.h)
Returns:
: a set of plane indicies satifying the request
Exceptions:
: geo::Exception if initialization fails

Definition at line 400 of file Geometry.cxx.

Referenced by mcchk::DetAna::Ana(), mcchk::CosmicAna::Ana(), vali::Validator::BookPT(), clusterss::MakeClusterSS::ClusterWindows(), clusterss::MakeClusterSS::FindLongWindow(), subshower::RecoSubShower::FindMaxWindow(), clusterss::MakeClusterSS::FindPlaneFromShwMax(), and testFindPlanes().

00401 {
00402   // Return the user's choice, fall back on all planes
00403   if      (v==kX) return fXplanes;
00404   else if (v==kY) return fYplanes;
00405   return fAllPlanes;
00406 }

const CellGeo & Geometry::IdToCell const CellUniqueId id,
int *  iplane = 0,
int *  icell = 0
const
 

Given a unique cell identifier, look up the plane and cell numbers.

Parameters:
id : Cell identifier
iplane : optional, pointer to return plane number
icell : optional, pointer to return cell number
Returns:
cell geometry description
Exceptions:
geo::Exception if id is not found

Definition at line 380 of file Geometry.cxx.

References geo::PlaneGeo::Cell(), fIdMap, and Plane().

Referenced by CurrentCell(), novamc::MCApplication::Stepping(), and testUniqueId().

00383 {
00384   UIDMap::const_iterator itr = fIdMap.find(id);
00385   if (itr == fIdMap.end()) {
00386     throw Exception(__FILE__,__LINE__,Exception::BAD_CELL_UNIQUE_ID);
00387   }
00388   *iplane = (int)itr->second.first;
00389   *icell  = (int)itr->second.second;
00390   return this->Plane(*iplane).Cell(*icell);
00391 }

Geometry & Geometry::Instance int  runnum,
const char *  gdmlfile
[static]
 

Build a detector geometry

Parameters:
runnum : Run number to configure geometry for
detId : Detector ID (0=near, 1=far, 2=ipnd)
Exceptions:
geo::Exception if requested detector/run combination cannot be satisfied

Definition at line 49 of file Geometry.cxx.

References fGDMLFile, fInstance, fRunNumber, and Geometry().

Referenced by clusterss::CompareClusters::Ana(), clusterss::ClusterCheck::Ana(), util::EDMUtils::GetGeometry(), main(), clusterss::MakeClusterSS::Reco(), evgen::CosmicsGen::Reco(), and cellhitmerge::CellHitMerge::Reco().

00050 {
00051   // If our current geometry doesn't match request, get rid of it
00052   if (fInstance!=0) {
00053     if (fInstance->fRunNumber != runnum || fInstance->fGDMLFile.compare(gdmlfile) != 0) {
00054       delete fInstance;
00055       fInstance = 0;
00056       gGeoManager->Clear();
00057     }
00058   }
00059 
00060   // Build a new geometry if we need to
00061   if (fInstance==0) {
00062     std::string gdml;
00063     if (gdmlfile!=0) {
00064       // User has specficied the file name to load
00065       gdml = gdmlfile;
00066     }
00067     else {
00068       char msg[256];
00069       sprintf(msg,"no geometry file specified");
00070       throw Exception(__FILE__,__LINE__,Exception::BAD_GEO_CONFIG,msg);
00071     }
00072     if (runnum<0) {
00073       char msg[256];
00074       std::sprintf(msg,"Bad run number: %d",runnum);
00075       throw Exception(__FILE__,__LINE__,Exception::BAD_GEO_CONFIG,msg);
00076     }
00077     fInstance = new Geometry(gdml.c_str());
00078     fInstance->fGDMLFile  = gdml;
00079     fInstance->fRunNumber = runnum;
00080   }//end if no instance
00081   return *fInstance;
00082 }

void Geometry::MakePlane std::vector< const TGeoNode * > &  n,
unsigned int  depth
[private]
 

Definition at line 303 of file Geometry.cxx.

References fPlanes.

Referenced by FindPlanes().

00305 {
00306   fPlanes.push_back(new PlaneGeo(n, depth));
00307 }

double Geometry::MassBetweenPoints double *  p1,
double *  p2
const
 

Return the column density between 2 points in gm/cm^2

Parameters:
p1 : pointer to array holding xyz of first point in world coordinates
p2 : pointer to array holding xyz of second point in world coorinates

Definition at line 531 of file Geometry.cxx.

00532 {
00533   // The purpose of this method is to determine the column density
00534   // between the two points given.  Do that by starting at p1 and 
00535   // stepping until you get to the node of p2.  calculate the distance
00536   // between the point just inside that node and p2 to get the last
00537   // bit of column density
00538   double columnD = 0.;
00539 
00540   // first initialize a track - get the direction cosines
00541   double length = TMath::Sqrt(TMath::Power(p2[0]-p1[0], 2.)
00542                               + TMath::Power(p2[1]-p1[1], 2.)
00543                               + TMath::Power(p2[2]-p1[2], 2.));
00544   double dxyz[3] = {(p2[0]-p1[0])/length, 
00545                     (p2[1]-p1[1])/length,
00546                     (p2[2]-p1[2])/length}; 
00547 
00548   gGeoManager->InitTrack(p1,dxyz);
00549   
00550   // might be helpful to have a point to a TGeoNode
00551   TGeoNode *node = gGeoManager->GetCurrentNode();
00552 
00553   // check that the points are not in the same volume already.  
00554   // if they are in different volumes, keep stepping until you 
00555   // are in the same volume as the second point
00556   while(!gGeoManager->IsSameLocation(p2[0], p2[1], p2[2])){
00557     gGeoManager->FindNextBoundary();
00558 
00559     columnD +=
00560       gGeoManager->GetStep() * 
00561       node->GetMedium()->GetMaterial()->GetDensity();
00562     
00563     // the act of stepping puts you in the next node and returns that
00564     // node
00565     node = gGeoManager->Step();
00566   } //end loop to get to volume of second point
00567 
00568   // now you are in the same volume as the last point, but not at that
00569   // point.  get the distance between the current point and the last
00570   // one
00571   const double *current = gGeoManager->GetCurrentPoint();
00572   length = TMath::Sqrt(TMath::Power(p2[0]-current[0], 2.)
00573                        + TMath::Power(p2[1]-current[1], 2.)
00574                        + TMath::Power(p2[2]-current[2], 2.));
00575   columnD += length*node->GetMedium()->GetMaterial()->GetDensity();
00576 
00577   return columnD;
00578 }

const TGeoMaterial * Geometry::Material double  x,
double  y,
double  z
const
 

Definition at line 182 of file Geometry.cxx.

00183 {
00184   const TGeoNode* node = gGeoManager->FindNode(x,y,z);
00185   return node->GetMedium()->GetMaterial();
00186 }

const unsigned int Geometry::NextPlaneInView unsigned int  p1,
int  d = +1
const
 

Return the index of the next plane in a particular view

Parameters:
p1 : Index of current plane
d : Direction to step in (>0 next down stream, <0 upstream)
Returns:
Index of next plane in same view. If no plane can be found (eg. reach end of detector) return kPLANE_NOT_FOUND.
Exceptions:
: Asserts that p1 is in valid range

Definition at line 419 of file Geometry.cxx.

References fPlanes.

Referenced by clusterss::MakeClusterSS::ClusterWindows(), subshower::RecoSubShower::LinkedBlobs(), subshower::RecoSubShower::NumXYPlanes(), calhit::CalHit::RoughReco(), clusterss::MakeClusterSS::StupidHoughNoWork(), testFindPlanes(), and subshower::RecoSubShower::TransCluster().

00420 {
00421   assert(p1<fPlanes.size());
00422 
00423   // March along from p1 in the direction d until we find a match
00424   int s = (d>0 ? 1 : -1);
00425   if (s<0 && p1==0)                return kPLANE_NOT_FOUND;
00426   if (s>0 && p1==fPlanes.size()-1) return kPLANE_NOT_FOUND;
00427   for (unsigned int p2=p1+d; 1; p2+=s) {
00428     if (fPlanes[p1]->View() == fPlanes[p2]->View()) return p2;
00429     if (p2==fPlanes.size()-1) return kPLANE_NOT_FOUND;
00430     if (p2==0)                return kPLANE_NOT_FOUND;
00431   }
00432 }

const unsigned int Geometry::NextPlaneOtherView unsigned int  p1,
int  d = +1
const
 

Return the index of the next plane in a particular view

Parameters:
p1 : Index of current plane
d : Direction to step in (>0 next down stream, <0 upstream)
Returns:
index of next plane in same view.
Exceptions:
: Asserts that p1 is in valid range

Definition at line 444 of file Geometry.cxx.

References fPlanes.

Referenced by calhit::CalHit::RoughReco(), and testFindPlanes().

00445 {
00446   assert(p1<fPlanes.size());
00447 
00448   // March along from p1 in the direction d until we find a match
00449   int s = (d>0 ? 1 : -1);
00450   if (s<0 && p1==0)                return kPLANE_NOT_FOUND;
00451   if (s>0 && p1==fPlanes.size()-1) return kPLANE_NOT_FOUND;
00452   for (unsigned int p2=p1+d; 1; p2+=s) {
00453     if (fPlanes[p1]->View() != fPlanes[p2]->View()) return p2;
00454     if (p2==fPlanes.size()) return kPLANE_NOT_FOUND;
00455     if (p2==0)              return kPLANE_NOT_FOUND;
00456   }
00457 }

unsigned int geo::Geometry::NPlanes  )  const [inline]
 

Definition at line 49 of file Geometry.h.

Referenced by CellInfo(), calib::Calibrator::GetXYZD(), and cellhitmerge::CellHitMerge::Reco().

00049 { return fPlanes.size(); }

const PlaneGeo & Geometry::Plane unsigned int  i  )  const
 

Return the geometry description of the ith plane in the detector.

Parameters:
i : input plane number, starting from 0
Returns:
plane geometry for ith plane
Exceptions:
geo::Exception if iplane is outside allowed range

Definition at line 318 of file Geometry.cxx.

References fPlanes.

Referenced by spider::Filament::AddHit(), mcchk::DetAna::Ana(), mcchk::CosmicAna::Ana(), subshower::RecoSubShower::BestHough(), clusterss::MakeClusterSS::BestHough(), vali::Validator::BookPT(), BuildMaps(), rpr::TrackReco::CalcChi2Ndof(), subshower::RecoSubShower::CalculateEnergyVertexAngle(), evd::RecoBaseDrawer::CellHit2D(), CellInfo(), clusterss::MakeClusterSS::ClusterWindows(), rpr::TrackReco::Connection(), evd::Display3D::DrawMCHits(), evd::Display3D::DrawRawData(), clusterss::MakeClusterSS::FindAllTransverseWindows(), spider::Filament::FitDist(), evd::SimulationDrawer::FLSHit2D(), calib::Calibrator::GetXYZD(), subshower::RecoSubShower::HoughTransCluster(), clusterss::MakeClusterSS::HoughTransCluster(), IdToCell(), main(), calib::Calibrator::MakeCellHit(), ctrk::CosmicTrack::MakeClusters(), ctrk::CosmicTrack::MakeTracks(), vali::Validator::PlotCell(), vali::Validator::PlotPT(), raw_digit_sort(), evd::RawDataDrawer::RawDigit2D(), rpr::TrackReco::Reco(), spider::SpiderWeb::Reco(), photrans::SimpleTransport::Reco(), rpr::FindTrackSeg::Reco(), trk::DemoShell::Reco(), cellhitmerge::CellHitMerge::Reco(), calhit::CalHit::RoughReco(), SetDetectorSize(), clusterss::MakeClusterSS::StupidHoughNoWork(), testCellPos(), testFindCell(), testFindPlanes(), testStepping(), testUniqueId(), subshower::RecoSubShower::TransCluster(), spider::Filament::Update(), and cluster::Clusterer::UpdateClusterStat().

00319 {
00320   if (i>=fPlanes.size()) throw Exception(__FILE__,__LINE__,
00321                                          Exception::BAD_PLANE_NUMBER,
00322                                          "Plane index out of range");
00323   return *fPlanes[i];
00324 }

TGeoManager * Geometry::ROOTGeoManager  )  const
 

Definition at line 135 of file Geometry.cxx.

Referenced by photrans::SimpleTransport::Reco().

00136 {
00137   return gGeoManager;
00138 }

void Geometry::SetDetectorSize  )  [private]
 

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 586 of file Geometry.cxx.

References geo::PlaneGeo::Cell(), fDetHalfHeight, fDetHalfWidth, fDetLength, fPlanes, geo::CellGeo::GetCenter(), and Plane().

Referenced by Geometry().

00587 {
00588   double dummy;
00589 
00590   // Compute the height and width based on the sizes of the vertical
00591   // and horizontal planes. Remember, in the local plane frame z goes
00592   // along the length of the modules. Hence the use of axis=3 to get
00593   // the detector height (from the vertical planes) and width (from
00594   // the horizontal planes).
00595   const TGeoVolume* vvp = gGeoManager->GetVolume("vPlaneV");
00596   const TGeoVolume* hvp = gGeoManager->GetVolume("vPlaneH");
00597 
00598   // Handle the case of the IPND where the volumes are names slightly
00599   // different accounting for the different widths of the modules.
00600   if (vvp==0) vvp = gGeoManager->GetVolume("vPlaneVIPND");
00601   if (hvp==0) hvp = gGeoManager->GetVolume("vPlaneHIPND");
00602   if (vvp==0 && hvp==0) {
00603     vvp = gGeoManager->GetVolume("vPlaneVIPND");
00604     hvp = gGeoManager->GetVolume("vPlaneHIPND");
00605   }
00606   if (vvp==0 || hvp==0) {
00607     throw Exception(__FILE__,__LINE__,
00608                     Exception::BAD_GEO_CONFIG,
00609                     "Unable to find shapes to set detector size");
00610   }
00611   
00612   const TGeoShape* vp = vvp->GetShape();
00613   const TGeoShape* hp = hvp->GetShape();
00614   vp->GetAxisRange(3,dummy,fDetHalfHeight);
00615   hp->GetAxisRange(3,dummy,fDetHalfWidth);
00616 
00617   // Remember -- in the local plane coordinate frame, the "depth" is
00618   // along y, hence axis=2. In the world coordinate this is along z.
00619   double vplanez1, vplanez2;
00620   vp->GetAxisRange(2,vplanez1,vplanez2);
00621   
00622   // Compute the detector length based on the center positions of
00623   // cells in the first and last plane. Adjust to account for the
00624   // plane widths. This ajustment works if the vertical and horizontal
00625   // planes have the same depth, or if the detector begins and ends
00626   // with vertical planes. As of the TDR both of these conditions was
00627   // true for all the NOvA detectors.
00628   double xyz1[3] = {0,0,0};
00629   double xyz2[3] = {0,0,0};
00630   this->Plane(0).               Cell(0).GetCenter(xyz1);
00631   this->Plane(fPlanes.size()-1).Cell(0).GetCenter(xyz2);
00632   fDetLength = (xyz2[2]-xyz1[2])+(vplanez2-vplanez1);
00633 }

void Geometry::SetDrawOptions  )  [private]
 

Definition at line 267 of file Geometry.cxx.

Referenced by Geometry().

00267 { }

double Geometry::SurfaceY  )  const
 

A typical y-position value at the surface (where earth meets air) for this detector site

Returns:
typical y position at surface in units of cm
Exceptions:
geo::Exception if detector ID is not set properly

Todo:
: Only far det number is reasonable accurate. Near and IPND need some thought

Definition at line 222 of file Geometry.cxx.

References fGDMLFile.

00223 {
00224   if(fGDMLFile.find("near")      != std::string::npos) return 130.0E2; // near: surface ~130m above det. center
00225   else if(fGDMLFile.find("far")  != std::string::npos) return 9.94E2;  // far : surface ~10m  above det. center
00226   else if(fGDMLFile.find("ipnd") != std::string::npos) return -1.966E2;// IPND: surface ~2m   below det. center
00227   else throw Exception(__FILE__,__LINE__,Exception::BAD_GEO_CONFIG);
00228  
00229   return 0.0;
00230 }

double Geometry::TotalMass  )  const
 

Return the total mass of the detector

Definition at line 496 of file Geometry.cxx.

Referenced by main().

00497 {
00498   double mass = 0.;
00499 
00500   // the detector resides in the vDetEnclosure volume so start there
00501   const TGeoVolume *enclosure = gGeoManager->FindVolumeFast("vDetEnclosure");
00502   
00503   //
00504   // loop over the daughters in the volume and let ROOT calculate the
00505   // mass in kg for you for the daughter nodes that have vBlock in
00506   // their names this only goes one level down, so no worries of
00507   // double counting non-superblock blocks
00508   //
00509   for(int d = 0; d < enclosure->GetNdaughters(); ++d){
00510     const TGeoNode* node = enclosure->GetNode(d);
00511     const char* nm = node->GetName();
00512     if(strncmp(nm,"vSuperBlock",11) == 0){
00513       std::string name = nm;
00514       name.erase(name.find('_'));
00515       mass += gGeoManager->FindVolumeFast(name.c_str())->Weight();
00516     }
00517   }
00518   return mass;
00519 }

void Geometry::WorldBox double *  xlo,
double *  xhi,
double *  ylo,
double *  yhi,
double *  zlo,
double *  zhi
const
 

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.

Parameters:
xlo : On return, lower bound on x positions
xhi : On return, upper bound on x positions
ylo : On return, lower bound on y positions
yhi : On return, upper bound on y positions
zlo : On return, lower bound on z positions
zhi : On return, upper bound on z positions

Definition at line 244 of file Geometry.cxx.

Referenced by testProject().

00247 {
00248   const TGeoShape* s = gGeoManager->GetVolume("vWorld")->GetShape();
00249   assert(s);
00250 
00251   if (xlo || xhi) {
00252     double x1, x2;
00253     s->GetAxisRange(1,x1,x2); if (xlo) *xlo = x1; if (xhi) *xhi = x2;
00254   }
00255   if (ylo || yhi) {
00256     double y1, y2;
00257     s->GetAxisRange(2,y1,y2); if (ylo) *ylo = y1; if (yhi) *yhi = y2;
00258   }
00259   if (zlo || zhi) {
00260     double z1, z2;
00261     s->GetAxisRange(3,z1,z2); if (zlo) *zlo = z1; if (zhi) *zhi = z2;
00262   }
00263 }


Member Data Documentation

std::set<unsigned int> geo::Geometry::fAllPlanes [private]
 

List of all planes.

Definition at line 86 of file Geometry.h.

Referenced by BuildMaps().

double geo::Geometry::fDetHalfHeight [private]
 

Detector 1/2 height (cm).

Definition at line 91 of file Geometry.h.

Referenced by SetDetectorSize().

double geo::Geometry::fDetHalfWidth [private]
 

Detector 1/2 width (cm).

Definition at line 92 of file Geometry.h.

Referenced by SetDetectorSize().

double geo::Geometry::fDetLength [private]
 

Detector 1/2 length (cm).

Definition at line 90 of file Geometry.h.

Referenced by SetDetectorSize().

std::string geo::Geometry::fGDMLFile [private]
 

gdml file holding the geometry

Definition at line 82 of file Geometry.h.

Referenced by Instance(), and SurfaceY().

UIDMap geo::Geometry::fIdMap [private]
 

Unique ID -> Plane,Cell.

Definition at line 84 of file Geometry.h.

Referenced by BuildMaps(), and IdToCell().

Geometry * Geometry::fInstance = 0 [static, private]
 

Instance of geometry.

Definition at line 28 of file Geometry.cxx.

Referenced by Instance().

std::vector<PlaneGeo*> geo::Geometry::fPlanes [private]
 

The detector planes.

Definition at line 83 of file Geometry.h.

Referenced by BuildMaps(), Geometry(), MakePlane(), NextPlaneInView(), NextPlaneOtherView(), Plane(), SetDetectorSize(), and ~Geometry().

int geo::Geometry::fRunNumber [private]
 

Run number of configuration.

Definition at line 81 of file Geometry.h.

Referenced by Instance().

std::set<unsigned int> geo::Geometry::fXplanes [private]
 

List of X measuring planes.

Definition at line 87 of file Geometry.h.

Referenced by BuildMaps().

std::set<unsigned int> geo::Geometry::fYplanes [private]
 

List of Y measuring planes.

Definition at line 88 of file Geometry.h.

Referenced by BuildMaps().


The documentation for this class was generated from the following files:
Generated on Sun Nov 22 04:45:31 2009 for NOvA Offline by  doxygen 1.3.9.1