From 521583d37d74a4d2a3bef8570de039ef66779e74 Mon Sep 17 00:00:00 2001 From: echenard Date: Tue, 28 Jul 2026 17:55:23 -0500 Subject: [PATCH 1/3] Refactoring of the caloriemter base code and bug fixes --- Analyses/src/CaloCalibAna_module.cc | 4 +- Analyses/src/CaloCalib_module.cc | 6 +- Analyses/src/DiskCal00_module.cc | 15 +- Analyses/src/ReadBack_module.cc | 2 +- CalPatRec/inc/CalHelixFinderAlg.hh | 10 +- CalPatRec/src/AgnosticHelixFinder_module.cc | 12 +- CalPatRec/src/CalHelixFinderAlg.cc | 8 +- CalPatRec/src/CalHelixFinder_module.cc | 4 +- CalPatRec/src/CalLineTimePeakFinder_module.cc | 6 +- CalPatRec/src/CalTimePeakFinder_module.cc | 4 +- CalPatRec/src/DeltaFinder_types.cc | 2 +- CaloCluster/src/CaloClusterFast_module.cc | 9 +- CaloCluster/src/CaloTrigger_module.cc | 12 +- CaloCluster/src/ClusterFinder.cc | 7 +- CaloCluster/src/ClusterUtils.cc | 8 +- CaloDiag/src/CaloExample_module.cc | 4 +- CaloDiag/src/CaloNNDiag_module.cc | 8 +- CaloDiag/src/CaloNeutron_module.cc | 10 +- CaloFilters/src/CaloCosmicCalib_module.cc | 4 +- CaloMC/inc/ShowerStepUtil.hh | 2 +- CaloMC/src/CaloDigiMaker_module.cc | 4 +- CaloMC/src/CaloPhotonPropagation.cc | 2 +- CaloMC/src/CaloShowerROMaker_module.cc | 2 +- CaloMC/src/CaloShowerStepMaker_module.cc | 30 +- CaloMC/src/ShowerStepUtil.cc | 2 +- CaloReco/src/CaloHitMaker_module.cc | 2 +- CalorimeterGeom/CMakeLists.txt | 10 +- CalorimeterGeom/data/crystalPos.txt | 1348 +++++++++++++++++ CalorimeterGeom/inc/CaloG4Info.hh | 65 + CalorimeterGeom/inc/CaloGeomUtil.hh | 63 - CalorimeterGeom/inc/CaloInfo.hh | 66 - CalorimeterGeom/inc/Calorimeter.hh | 68 +- CalorimeterGeom/inc/Crystal.hh | 105 +- CalorimeterGeom/inc/CrystalMapper.hh | 41 +- CalorimeterGeom/inc/Disk.hh | 100 +- CalorimeterGeom/inc/DiskCalorimeter.hh | 81 +- CalorimeterGeom/inc/DiskGeomInfo.hh | 80 - CalorimeterGeom/inc/DiskInfo.hh | 106 ++ CalorimeterGeom/inc/SquareMapper.hh | 108 +- CalorimeterGeom/inc/SquareShiftMapper.hh | 107 +- CalorimeterGeom/src/CaloGeomUtil.cc | 122 -- CalorimeterGeom/src/Disk.cc | 585 +++---- CalorimeterGeom/src/DiskCalorimeter.cc | 202 +-- CalorimeterGeom/src/DiskInfo.cc | 54 + CalorimeterGeom/src/SConscript | 3 + CalorimeterGeom/src/SquareMapper.cc | 275 ++-- CalorimeterGeom/src/SquareShiftMapper.cc | 320 ++-- EventDisplay/src/DataInterface.cc | 36 +- .../src/CORSIKAEventGenerator_module.cc | 8 +- EventGenerator/src/CaloCalibGun_module.cc | 28 +- EventGenerator/src/CaloTBGun_module.cc | 2 +- EventGenerator/src/CosmicCRY.cc | 2 +- Filters/src/CaloDtsClusterFilter_module.cc | 6 +- Filters/src/CaloPileupDtsFilter_module.cc | 4 +- Filters/src/PhotonFilter_module.cc | 12 +- GeometryService/src/DiskCalorimeterMaker.cc | 434 +++--- GeometryService/src/VirtualDetectorMaker.cc | 15 +- Mu2eG4/geom/calorimeter_CsI_v2.txt | 1 + Mu2eG4/src/CaloReadoutCardSD.cc | 2 +- Mu2eG4/src/CaloReadoutSD.cc | 2 +- Mu2eG4/src/constructDiskCalorimeter.cc | 230 +-- Mu2eKinKal/inc/KKFit.hh | 14 +- .../src/RobustMultiHelixFinder_module.cc | 12 +- TrkReco/src/RobustHelixFit.cc | 2 +- TrkReco/src/TrkTimeCalculator.cc | 2 +- 65 files changed, 3015 insertions(+), 1885 deletions(-) create mode 100644 CalorimeterGeom/data/crystalPos.txt create mode 100644 CalorimeterGeom/inc/CaloG4Info.hh delete mode 100644 CalorimeterGeom/inc/CaloGeomUtil.hh delete mode 100644 CalorimeterGeom/inc/CaloInfo.hh delete mode 100644 CalorimeterGeom/inc/DiskGeomInfo.hh create mode 100644 CalorimeterGeom/inc/DiskInfo.hh delete mode 100644 CalorimeterGeom/src/CaloGeomUtil.cc create mode 100644 CalorimeterGeom/src/DiskInfo.cc diff --git a/Analyses/src/CaloCalibAna_module.cc b/Analyses/src/CaloCalibAna_module.cc index 85f6f86edd..a34dc235ab 100644 --- a/Analyses/src/CaloCalibAna_module.cc +++ b/Analyses/src/CaloCalibAna_module.cc @@ -304,7 +304,7 @@ namespace mu2e { const CaloHit& hit = CaloHits.at(ic); int diskId = cal.crystal(hit.crystalID()).diskID(); - CLHEP::Hep3Vector crystalPos = cal.geomUtil().mu2eToDiskFF(diskId,cal.crystal(hit.crystalID()).position()); //in disk FF frame + CLHEP::Hep3Vector crystalPos = cal.mu2eToDiskFF(diskId,cal.crystal(hit.crystalID()).position()); //in disk FF frame const auto eDepMCs = caloHitTruth[ic].energyDeposits(); @@ -492,7 +492,7 @@ namespace mu2e { if ( (hit.volumeId()VirtualDetectorId::EMC_Disk_1_EdgeOut) && hit.volumeId() != VirtualDetectorId::TT_Back) continue; - CLHEP::Hep3Vector VDPos = cal.geomUtil().mu2eToTracker(hit.position()); + CLHEP::Hep3Vector VDPos = cal.mu2eToTracker(hit.position()); vdId_[nVd_] = hit.volumeId(); vdPdgId_[nVd_] = hit.simParticle()->pdgId(); diff --git a/Analyses/src/CaloCalib_module.cc b/Analyses/src/CaloCalib_module.cc index e8c8c36b62..a0bb21cba4 100644 --- a/Analyses/src/CaloCalib_module.cc +++ b/Analyses/src/CaloCalib_module.cc @@ -253,7 +253,7 @@ namespace mu2e { //CaloHitCollection const& CaloHits(*CaloHitsHandle); //const double CrDensity = 4.9e-6; // in kg/mm3 to be consistent with volume units! - //const double CrMass = CrDensity*cal.caloInfo().crystalVolume(); + //const double CrMass = CrDensity*cal.G4Info().crystalVolume(); @@ -290,7 +290,7 @@ namespace mu2e { { CaloHit const& hit = CaloHits.at(ic); int diskId = cal.crystal(hit.id()).diskID(); - CLHEP::Hep3Vector crystalPos = cal.geomUtil().mu2eToDiskFF(diskId,cal.crystal(hit.id()).position()); //in disk FF frame + CLHEP::Hep3Vector crystalPos = cal.mu2eToDiskFF(diskId,cal.crystal(hit.id()).position()); //in disk FF frame CaloHit const& caloHit = *(hit.readouts().at(0)); @@ -319,7 +319,7 @@ namespace mu2e { while (grandMother->hasParent()) grandMother = grandMother->parent(); GenParticle const* generated = grandMother->genParticle() ? grandMother->genParticle().get() : 0; - CLHEP::Hep3Vector hitSimPos = cal.geomUtil().mu2eToDiskFF(diskId,hitSim.position().at(ip)); //in disk FF frame + CLHEP::Hep3Vector hitSimPos = cal.mu2eToDiskFF(diskId,hitSim.position().at(ip)); //in disk FF frame _motId[_nSim] = mother->id().asInt(); _motPdgId[_nSim] = mother->pdgId(); diff --git a/Analyses/src/DiskCal00_module.cc b/Analyses/src/DiskCal00_module.cc index 78ac4ea2df..f2aa4fc8a6 100644 --- a/Analyses/src/DiskCal00_module.cc +++ b/Analyses/src/DiskCal00_module.cc @@ -19,7 +19,7 @@ #include "art_root_io/TFileService.h" // Mu2e includes. -#include "Offline/CalorimeterGeom/inc/DiskCalorimeter.hh" +#include "Offline/CalorimeterGeom/inc/Calorimeter.hh" #include "Offline/DataProducts/inc/CaloConst.hh" #include "Offline/GeometryService/inc/GeomHandle.hh" #include "Offline/RecoDataProducts/inc/CaloHit.hh" @@ -155,7 +155,7 @@ namespace mu2e { _hNHits->Fill(cryHits.size()); maxHits = std::max( cryHits.size(), maxHits ); - GeomHandle calGeom; + GeomHandle calGeom; double totalEdep(0.); set hit_crystals; @@ -265,18 +265,17 @@ namespace mu2e { } void DiskCal00::printCalInfo(){ - DiskCalorimeter const& cal(*GeomHandle()); - int nSiPM = cal.nCrystals()*cal.caloInfo().getInt("nSiPMPerCrystal"); + Calorimeter const& cal(*GeomHandle()); + int nSiPM = cal.nCrystals()*cal.G4Info().get("nSiPMPerCrystal"); cout << "Information about the disk Calorimeter: " << endl; cout << "Number of disks: " << cal.nDisks() << endl; cout << "Number of Readouts: " << nSiPM << " " << CaloConst::_nSiPMPerCrystal << " " << nSiPM/CaloConst::_nSiPMPerCrystal << endl; - cout << "Hex side size: " << 2.0*cal.caloInfo().getDouble("crystalXYLength") << endl; + cout << "Hex side size: " << 2.0*cal.G4Info().get("crystalXYLength") << endl; - cout << "Depth: " << cal.caloInfo().getDouble("crystalZLength") << endl; - cout << "Origin: " << cal.geomUtil().origin() << endl; + cout << "Depth: " << cal.G4Info().get("crystalZLength") << endl; for (unsigned i=0; iFill(hitPos.x(),hitPos.y()); if ( _diagLevel > 1 && _nAnalyzed < _maxFullPrint ) diff --git a/CalPatRec/inc/CalHelixFinderAlg.hh b/CalPatRec/inc/CalHelixFinderAlg.hh index 1fdf230fad..b47d7c2516 100644 --- a/CalPatRec/inc/CalHelixFinderAlg.hh +++ b/CalPatRec/inc/CalHelixFinderAlg.hh @@ -327,11 +327,11 @@ namespace mu2e { void setTracker (const Tracker* Tracker) { _tracker = Tracker; } void setCalorimeter(const Calorimeter* Cal ) { _calorimeter = Cal; - fCaloOffset = (_calorimeter->caloInfo().getDouble("diskCaseZLength")/2. - + ( _calorimeter->caloInfo().getDouble("BPPipeZOffset") - + _calorimeter->caloInfo().getDouble("BPHoleZLength") - + _calorimeter->caloInfo().getDouble("FEEZLength"))/2. - - _calorimeter->caloInfo().getDouble("FPCarbonZLength") - _calorimeter->caloInfo().getDouble("FPFoamZLength") + fCaloOffset = (_calorimeter->G4Info().get("diskCaseZLength")/2. + + ( _calorimeter->G4Info().get("BPPipeZOffset") + + _calorimeter->G4Info().get("BPHoleZLength") + + _calorimeter->G4Info().get("FEEZLength"))/2. + - _calorimeter->G4Info().get("FPCarbonZLength") - _calorimeter->G4Info().get("FPFoamZLength") ); } //----------------------------------------------------------------------------- diff --git a/CalPatRec/src/AgnosticHelixFinder_module.cc b/CalPatRec/src/AgnosticHelixFinder_module.cc index 7135d8745f..ef711554c7 100644 --- a/CalPatRec/src/AgnosticHelixFinder_module.cc +++ b/CalPatRec/src/AgnosticHelixFinder_module.cc @@ -405,10 +405,10 @@ namespace mu2e { _bz0ConvSq = _bz0Conv * _bz0Conv; // Offset for calo cluster z positions - double offset = _calorimeter->caloInfo().getDouble("diskCaseZLength"); - offset += _calorimeter->caloInfo().getDouble("BPPipeZOffset"); - offset += _calorimeter->caloInfo().getDouble("BPHoleZLength"); - offset += _calorimeter->caloInfo().getDouble("FEEZLength"); + double offset = _calorimeter->G4Info().get("diskCaseZLength"); + offset += _calorimeter->G4Info().get("BPPipeZOffset"); + offset += _calorimeter->G4Info().get("BPHoleZLength"); + offset += _calorimeter->G4Info().get("FEEZLength"); offset /= 2.0; _caloZOffset = offset; } @@ -680,8 +680,8 @@ namespace mu2e { if (cl.isNonnull()) { cHit hit; hit.hitIndice = HitType::CALOCLUSTER; - CLHEP::Hep3Vector gpos = _calorimeter->geomUtil().diskToMu2e(cl->diskID(), cl->cog3Vector()); - CLHEP::Hep3Vector tpos = _calorimeter->geomUtil().mu2eToTracker(gpos); + CLHEP::Hep3Vector gpos = _calorimeter->diskToMu2e(cl->diskID(), cl->cog3Vector()); + CLHEP::Hep3Vector tpos = _calorimeter->mu2eToTracker(gpos); _caloPos.SetCoordinates(tpos.x(), tpos.y(), tpos.z() - _caloZOffset); hit.pos = _caloPos; _tcHits.push_back(hit); diff --git a/CalPatRec/src/CalHelixFinderAlg.cc b/CalPatRec/src/CalHelixFinderAlg.cc index b0c8d3400e..10749820d8 100644 --- a/CalPatRec/src/CalHelixFinderAlg.cc +++ b/CalPatRec/src/CalHelixFinderAlg.cc @@ -252,8 +252,8 @@ namespace mu2e { return; } //fill the calorimeter cluster info - Hep3Vector gpos = _calorimeter->geomUtil().diskToMu2e(cl->diskID(),cl->cog3Vector()); - Hep3Vector tpos = _calorimeter->geomUtil().mu2eToTracker(gpos); + Hep3Vector gpos = _calorimeter->diskToMu2e(cl->diskID(),cl->cog3Vector()); + Hep3Vector tpos = _calorimeter->mu2eToTracker(gpos); fCaloTime = cl->time(); fCaloX = tpos.x(); fCaloY = tpos.y(); @@ -272,8 +272,8 @@ namespace mu2e { // if (cl == NULL) return false; // //fill the calorimeter cluster info - // Hep3Vector gpos = _calorimeter->geomUtil().diskToMu2e(cl->diskID(),cl->cog3Vector()); - // Hep3Vector tpos = _calorimeter->geomUtil().mu2eToTracker(gpos); + // Hep3Vector gpos = _calorimeter->diskToMu2e(cl->diskID(),cl->cog3Vector()); + // Hep3Vector tpos = _calorimeter->mu2eToTracker(gpos); // fCaloTime = cl->time(); // fCaloX = tpos.x(); // fCaloY = tpos.y(); diff --git a/CalPatRec/src/CalHelixFinder_module.cc b/CalPatRec/src/CalHelixFinder_module.cc index ecad489323..f0064cc070 100644 --- a/CalPatRec/src/CalHelixFinder_module.cc +++ b/CalPatRec/src/CalHelixFinder_module.cc @@ -416,9 +416,9 @@ namespace mu2e { if(HfResult._timeClusterPtr->hasCaloCluster()) { CLHEP::Hep3Vector gpos,tpos; - gpos = _hfinder._calorimeter->geomUtil().diskToMu2e(HfResult._timeClusterPtr->caloCluster()->diskID(), + gpos = _hfinder._calorimeter->diskToMu2e(HfResult._timeClusterPtr->caloCluster()->diskID(), HfResult._timeClusterPtr->caloCluster()->cog3Vector()); - tpos = _hfinder._calorimeter->geomUtil().mu2eToTracker(gpos); + tpos = _hfinder._calorimeter->mu2eToTracker(gpos); hel_t0 = HfResult._timeClusterPtr->caloCluster()->time() - (tpos.z() - z0)/sin(pitchAngle)/(beta*CLHEP::c_light); diff --git a/CalPatRec/src/CalLineTimePeakFinder_module.cc b/CalPatRec/src/CalLineTimePeakFinder_module.cc index 8322777b93..af80b2308f 100644 --- a/CalPatRec/src/CalLineTimePeakFinder_module.cc +++ b/CalPatRec/src/CalLineTimePeakFinder_module.cc @@ -135,11 +135,11 @@ namespace mu2e { stopping_target_ = sth.get(); // get the offset between the calo disks and the tracker system - calo_d0_offset_ = calorimeter_->geomUtil().mu2eToTracker(calorimeter_->geomUtil().diskToMu2e(0, CLHEP::Hep3Vector(0., 0., 0.))).z(); - calo_d1_offset_ = calorimeter_->geomUtil().mu2eToTracker(calorimeter_->geomUtil().diskToMu2e(1, CLHEP::Hep3Vector(0., 0., 0.))).z(); + calo_d0_offset_ = calorimeter_->mu2eToTracker(calorimeter_->diskToMu2e(0, CLHEP::Hep3Vector(0., 0., 0.))).z(); + calo_d1_offset_ = calorimeter_->mu2eToTracker(calorimeter_->diskToMu2e(1, CLHEP::Hep3Vector(0., 0., 0.))).z(); // get the target position in the tracker system - target_pos_ = calorimeter_->geomUtil().mu2eToTracker(stopping_target_->centerInMu2e()); + target_pos_ = calorimeter_->mu2eToTracker(stopping_target_->centerInMu2e()); } diff --git a/CalPatRec/src/CalTimePeakFinder_module.cc b/CalPatRec/src/CalTimePeakFinder_module.cc index 9b88f717c9..81a331065b 100644 --- a/CalPatRec/src/CalTimePeakFinder_module.cc +++ b/CalPatRec/src/CalTimePeakFinder_module.cc @@ -200,8 +200,8 @@ namespace mu2e { // convert cluster coordinates defined in the disk frame to the detector // coordinate system //----------------------------------------------------------------------------- - gpos = _calorimeter->geomUtil().diskToMu2e(cl->diskID(),cl->cog3Vector()); - tpos = _calorimeter->geomUtil().mu2eToTracker(gpos); + gpos = _calorimeter->diskToMu2e(cl->diskID(),cl->cog3Vector()); + tpos = _calorimeter->mu2eToTracker(gpos); xcl = tpos.x(); ycl = tpos.y(); diff --git a/CalPatRec/src/DeltaFinder_types.cc b/CalPatRec/src/DeltaFinder_types.cc index d0b5475347..a93e86e171 100644 --- a/CalPatRec/src/DeltaFinder_types.cc +++ b/CalPatRec/src/DeltaFinder_types.cc @@ -94,7 +94,7 @@ namespace mu2e { mu2e::GeomHandle tH; tracker = tH.get(); - // mu2e::GeomHandle cH; + // mu2e::GeomHandle cH; // calorimeter = cH.get(); //----------------------------------------------------------------------------- diff --git a/CaloCluster/src/CaloClusterFast_module.cc b/CaloCluster/src/CaloClusterFast_module.cc index 11496fad55..baf40613d5 100644 --- a/CaloCluster/src/CaloClusterFast_module.cc +++ b/CaloCluster/src/CaloClusterFast_module.cc @@ -118,7 +118,7 @@ namespace mu2e { ++iterStart; //start the clustering algorithm for the hits between iStart and iStop - std::queue crystalToVisit; + std::queue crystalToVisit; std::vector isVisited(cal.nCrystals()); //put the first hit in the cluster list @@ -134,7 +134,8 @@ namespace mu2e { isVisited[visitId]=true; auto neighborsId = cal.crystal(visitId).neighbors(); - if (extendSearch_) std::copy(cal.nextNeighbors(visitId).begin(), cal.nextNeighbors(visitId).end(), std::back_inserter(neighborsId)); + auto nn = cal.crystal(visitId).nextNeighbors(); + if (extendSearch_) std::copy(nn.begin(), nn.end(), std::back_inserter(neighborsId)); for (const auto& iId : neighborsId) { @@ -145,7 +146,7 @@ namespace mu2e { for (auto it=iterStart; it != iterStop; ++it) { if (*it==hits.size()) continue; - if (caloHits[*it].crystalID() != iId) continue; + if (caloHits[*it].crystalID() != static_cast(iId)) continue; if (caloHits[*it].energyDep() > ExpandCut_) crystalToVisit.push(iId); clusterList.push_back(*it); @@ -155,7 +156,7 @@ namespace mu2e { crystalToVisit.pop(); } - auto functorEnergy = [&caloHits](int a, int b) {return caloHits[a].energyDep() > caloHits[b].energyDep();}; + auto functorEnergy = [&caloHits](auto a, auto b) {return caloHits[a].energyDep() > caloHits[b].energyDep();}; std::sort(clusterList.begin(),clusterList.end(),functorEnergy); fillCluster(cal, caloHitsHandle, caloHits, clusterList, caloClusters ); diff --git a/CaloCluster/src/CaloTrigger_module.cc b/CaloCluster/src/CaloTrigger_module.cc index 33106486fa..ab67e5ecfe 100644 --- a/CaloCluster/src/CaloTrigger_module.cc +++ b/CaloCluster/src/CaloTrigger_module.cc @@ -134,7 +134,7 @@ namespace mu2e { mu2e::GeomHandle ch; const Calorimeter* cal = ch.get(); - int nro = cal->caloInfo().getInt("nSiPMPerCrystal"); + int nro = cal->G4Info().get("nSiPMPerCrystal"); unsigned offsetT0_ = unsigned(blindTime_/digiSampling_); unsigned nBinTime = unsigned (mbtime_ - blindTime_ + endTimeBuffer_) / digiSampling_; @@ -192,7 +192,7 @@ namespace mu2e { int ring1max(0); int ring1max2(0); std::queue crystalRing1; - for (const auto& nid : cal->neighbors(seed->crId_)) crystalRing1.push(nid); + for (const auto& nid : cal->crystal(seed->crId_).neighbors()) crystalRing1.push(nid); while (!crystalRing1.empty()){ int nid = crystalRing1.front(); for (int itimebin= seed->index_-DNTBINs_;itimebin<=seed->index_+DNTBINs_;++itimebin){ @@ -214,7 +214,7 @@ namespace mu2e { // int ring2max(0); std::queue crystalRing2; - for (const auto& nid : cal->nextNeighbors(seed->crId_)) crystalRing2.push(nid); + for (const auto& nid : cal->crystal(seed->crId_).nextNeighbors()) crystalRing2.push(nid); while (!crystalRing2.empty()){ int nid = crystalRing2.front(); for (int itimebin= seed->index_-DNTBINs_;itimebin<=seed->index_+DNTBINs_;++itimebin){ @@ -233,7 +233,7 @@ namespace mu2e { std::queue crystalToVisit; crystalToVisit.push(seed->crId_); - for (const auto& nid : cal->neighbors(seed->crId_)) crystalToVisit.push(nid); + for (const auto& nid : cal->crystal(seed->crId_).neighbors()) crystalToVisit.push(nid); if (diagLevel_ > 1) std::cout << " To visit: " << crystalToVisit.size() << std::endl; while (!crystalToVisit.empty()){ @@ -248,8 +248,8 @@ namespace mu2e { yc += cal->crystal(nid).localPosition().y()*hit.val_; hit.val_=-abs(hit.val_); if (diagLevel_ > 1) std::cout << " index add close to: " << nid << std::endl; - for (const auto& neighbor : cal->neighbors(nid)) crystalToVisit.push(neighbor); - if (extendSecond_) for (const auto& nneighbor : cal->nextNeighbors(nid)) crystalToVisit.push(nneighbor); + for (const auto& neighbor : cal->crystal(nid).neighbors()) crystalToVisit.push(neighbor); + if (extendSecond_) for (const auto& nneighbor : cal->crystal(nid).nextNeighbors()) crystalToVisit.push(nneighbor); } } crystalToVisit.pop(); diff --git a/CaloCluster/src/ClusterFinder.cc b/CaloCluster/src/ClusterFinder.cc index 6edad665ff..28dd0ee4b8 100644 --- a/CaloCluster/src/ClusterFinder.cc +++ b/CaloCluster/src/ClusterFinder.cc @@ -29,11 +29,12 @@ namespace mu2e { while (!crystalToVisit_.empty()) { - int visitId = crystalToVisit_.front(); + auto visitId = crystalToVisit_.front(); isVisited_[visitId] = true; - std::vector neighborsId = cal_->crystal(visitId).neighbors(); - if (addSecondRing_) neighborsId.insert(neighborsId.end(), cal_->nextNeighbors(visitId).begin(), cal_->nextNeighbors(visitId).end()); + auto neighborsId = cal_->crystal(visitId).neighbors(); + const auto& nn = cal_->crystal(visitId).nextNeighbors(); + if (addSecondRing_) neighborsId.insert(neighborsId.end(), nn.begin(), nn.end()); for (auto& iId : neighborsId) { diff --git a/CaloCluster/src/ClusterUtils.cc b/CaloCluster/src/ClusterUtils.cc index aecbae1b77..8ac3112ef3 100644 --- a/CaloCluster/src/ClusterUtils.cc +++ b/CaloCluster/src/ClusterUtils.cc @@ -49,9 +49,9 @@ namespace mu2e { double sx(0),sy(0),sz(0),sx2(0),sy2(0),sw(0); fill(sx,sy,sz,sx2,sy2,sw); - int iSection = cal_.crystal(hits_[0]->crystalID()).diskID(); + auto iSection = cal_.crystal(hits_[0]->crystalID()).diskID(); CLHEP::Hep3Vector cogMu2eFrame(sx/sw,sy/sw,sz/sw); - return cal_.geomUtil().mu2eToDiskFF(iSection,cogMu2eFrame); + return cal_.mu2eToDiskFF(iSection,cogMu2eFrame); } //------------------------------------------------------------------------------------------------- @@ -65,10 +65,10 @@ namespace mu2e { //------------------------------------------------------------------------------------------------- void ClusterUtils::fill(double& sx, double& sy, double& sz, double& sx2, double& sy2, double& sw) const { - int iSection = cal_.crystal(hits_[0]->crystalID()).diskID(); + auto iSection = cal_.crystal(hits_[0]->crystalID()).diskID(); for (const auto& hit : hits_) { - int crId(hit->crystalID()); + auto crId(hit->crystalID()); double energy(hit->energyDep()); if (cal_.crystal(crId).diskID() != iSection) continue; diff --git a/CaloDiag/src/CaloExample_module.cc b/CaloDiag/src/CaloExample_module.cc index 9e38a7329d..a743d5d0a8 100644 --- a/CaloDiag/src/CaloExample_module.cc +++ b/CaloDiag/src/CaloExample_module.cc @@ -282,7 +282,7 @@ namespace mu2e { { const CaloHit& hit = CaloHits.at(ic); int diskId = cal.crystal(hit.crystalID()).diskID(); - CLHEP::Hep3Vector crystalPos = cal.geomUtil().mu2eToDiskFF(diskId,cal.crystal(hit.crystalID()).position()); //in disk FF frame + CLHEP::Hep3Vector crystalPos = cal.mu2eToDiskFF(diskId,cal.crystal(hit.crystalID()).position()); //in disk FF frame const auto eDepMCs = caloHitTruth[ic].energyDeposits(); @@ -455,7 +455,7 @@ namespace mu2e { if ( (hit.volumeId()VirtualDetectorId::EMC_Disk_1_EdgeOut) && hit.volumeId() != VirtualDetectorId::TT_Back) continue; - CLHEP::Hep3Vector VDPos = cal.geomUtil().mu2eToTracker(hit.position()); + CLHEP::Hep3Vector VDPos = cal.mu2eToTracker(hit.position()); vdId_[nVd_] = hit.volumeId(); vdPdgId_[nVd_] = hit.simParticle()->pdgId(); diff --git a/CaloDiag/src/CaloNNDiag_module.cc b/CaloDiag/src/CaloNNDiag_module.cc index e97c366d71..c8442f40dc 100644 --- a/CaloDiag/src/CaloNNDiag_module.cc +++ b/CaloDiag/src/CaloNNDiag_module.cc @@ -280,7 +280,7 @@ namespace mu2e { for (unsigned int ic=0; icsecond->momentum().mag(); - simPos = cal.geomUtil().mu2eToDiskFF(cluster.diskID(), vdMapEntry->second->position()); + simPos = cal.mu2eToDiskFF(cluster.diskID(), vdMapEntry->second->position()); } sId.push_back(sim->id().asInt()); @@ -480,7 +480,7 @@ namespace mu2e { double hitTime = fmod(hit.time(),_mbtime); - CLHEP::Hep3Vector VDPos = cal.geomUtil().mu2eToTracker(hit.position()); + CLHEP::Hep3Vector VDPos = cal.mu2eToTracker(hit.position()); vdId_[nVd_] = hit.volumeId(); vdPdgId_[nVd_] = hit.simParticle()->pdgId(); diff --git a/CaloDiag/src/CaloNeutron_module.cc b/CaloDiag/src/CaloNeutron_module.cc index 52f49c3cbb..9488eca545 100644 --- a/CaloDiag/src/CaloNeutron_module.cc +++ b/CaloDiag/src/CaloNeutron_module.cc @@ -278,7 +278,7 @@ namespace mu2e { for (const auto& step : steps) { - CLHEP::Hep3Vector Pos = cal.geomUtil().mu2eToTracker(step.position()); + CLHEP::Hep3Vector Pos = cal.mu2eToTracker(step.position()); _cryId[_nCry] = step.volumeId(); _cryEdep[_nCry] = step.totalEDep(); @@ -299,7 +299,7 @@ namespace mu2e { for (const auto& step : steps) { - CLHEP::Hep3Vector Pos = cal.geomUtil().mu2eToTracker(step.position()); + CLHEP::Hep3Vector Pos = cal.mu2eToTracker(step.position()); _ROId[_nRO]=step.volumeId(); _ROEdep[_nRO]=step.totalEDep(); _ROPosX[_nRO]=Pos.x(); @@ -317,7 +317,7 @@ namespace mu2e { for (const auto& step : steps) { - CLHEP::Hep3Vector Pos = cal.geomUtil().mu2eToTracker(step.position()); + CLHEP::Hep3Vector Pos = cal.mu2eToTracker(step.position()); _ROCardId[_nROCard]=step.volumeId(); _ROCardEdep[_nROCard]=step.totalEDep(); _ROCardPosX[_nROCard]=Pos.x(); @@ -335,7 +335,7 @@ namespace mu2e { for (const auto& step : steps) { - CLHEP::Hep3Vector Pos = cal.geomUtil().mu2eToTracker(step.position()); + CLHEP::Hep3Vector Pos = cal.mu2eToTracker(step.position()); _CrateId[_nCrate]=step.volumeId(); _CrateEdep[_nCrate]=step.totalEDep(); _CratePosX[_nCrate]=Pos.x(); @@ -363,7 +363,7 @@ namespace mu2e { double hitTime = fmod(hit.time(),_mbtime); - CLHEP::Hep3Vector VDPos = cal.geomUtil().mu2eToTracker(hit.position()); + CLHEP::Hep3Vector VDPos = cal.mu2eToTracker(hit.position()); //CLHEP::Hep3Vector VDPos = cal.toTrackerFrame(hit.position()); _vdId[_nVd] = hit.volumeId(); diff --git a/CaloFilters/src/CaloCosmicCalib_module.cc b/CaloFilters/src/CaloCosmicCalib_module.cc index 4e9705cf9a..c0de796deb 100644 --- a/CaloFilters/src/CaloCosmicCalib_module.cc +++ b/CaloFilters/src/CaloCosmicCalib_module.cc @@ -100,7 +100,7 @@ namespace mu2e bool CaloCosmicCalib::beginRun(art::Run & run){ GeomHandle ch; _calogeom = ch.get(); - _crystalsize = _calogeom->caloInfo().getDouble("crystalXYLength") + 2.*_calogeom->caloInfo().getDouble("wrapperThickness"); + _crystalsize = _calogeom->G4Info().get("crystalXYLength") + 2.*_calogeom->G4Info().get("wrapperThickness"); return true; } @@ -152,7 +152,7 @@ namespace mu2e int _cryid = hitcalo->crystalID(); int DiskId = _calogeom->crystal(_cryid).diskID(); - CLHEP::Hep3Vector crystalpos = _calogeom->geomUtil().mu2eToDiskFF(DiskId,_calogeom->crystal(_cryid).position()); + CLHEP::Hep3Vector crystalpos = _calogeom->mu2eToDiskFF(DiskId,_calogeom->crystal(_cryid).position()); if( ene> _mincelcut){ // Use cells with energy above a MIP-like threshold xval = crystalpos.x(); diff --git a/CaloMC/inc/ShowerStepUtil.hh b/CaloMC/inc/ShowerStepUtil.hh index a080fd0b83..23a23f0555 100644 --- a/CaloMC/inc/ShowerStepUtil.hh +++ b/CaloMC/inc/ShowerStepUtil.hh @@ -19,7 +19,7 @@ namespace mu2e { w_(imax,0),pos_(0,0,0) {}; - void add(unsigned i, double eDepG4, double eDepVis, double time, double momentum, CLHEP::Hep3Vector& pos); + void add(unsigned i, double eDepG4, double eDepVis, double time, double momentum, const CLHEP::Hep3Vector& pos); void reset(unsigned i); void printBucket(unsigned i); diff --git a/CaloMC/src/CaloDigiMaker_module.cc b/CaloMC/src/CaloDigiMaker_module.cc index 8ec697cee0..dac577425c 100644 --- a/CaloMC/src/CaloDigiMaker_module.cc +++ b/CaloMC/src/CaloDigiMaker_module.cc @@ -192,14 +192,14 @@ namespace mu2e { mu2e::GeomHandle ch; calorimeter_ = ch.get(); - if (calorimeter_->nCrystals()<1 || calorimeter_->caloInfo().getInt("nSiPMPerCrystal")<1) return; + if (calorimeter_->nCrystals()<1 || calorimeter_->G4Info().get("nSiPMPerCrystal")<1) return; int waveformSize = (digitizationEnd_ - digitizationStart_ + startTimeBuffer_) / digiSampling_; if (ewMarker.spillType() != EventWindowMarker::SpillType::onspill) { waveformSize = (ewMarker.eventLength() - digitizationStart_ + startTimeBuffer_) / digiSampling_; } - int nWaveforms = calorimeter_->nCrystals()*calorimeter_->caloInfo().getInt("nSiPMPerCrystal"); + int nWaveforms = calorimeter_->nCrystals()*calorimeter_->G4Info().get("nSiPMPerCrystal"); if (waveformSize<1) throw cet::exception("Rethrow")<< "[CaloMC/CaloDigiMaker] digitization size too short " << std::endl; bool resetWaveform(false); std::vector waveform(waveformSize,0.0); diff --git a/CaloMC/src/CaloPhotonPropagation.cc b/CaloMC/src/CaloPhotonPropagation.cc index 7f28472a0e..d0165a128f 100644 --- a/CaloMC/src/CaloPhotonPropagation.cc +++ b/CaloMC/src/CaloPhotonPropagation.cc @@ -59,7 +59,7 @@ namespace mu2e { } const Calorimeter& cal = *(GeomHandle()); - lightSpeed_ = 300.0 / cal.caloInfo().getDouble("refractiveIndex"); //in mm/ns + lightSpeed_ = 300.0 / cal.G4Info().get("refractiveIndex"); //in mm/ns } //---------------------------------------------------------------------------- diff --git a/CaloMC/src/CaloShowerROMaker_module.cc b/CaloMC/src/CaloShowerROMaker_module.cc index 9dddc4a780..29b48d4138 100644 --- a/CaloMC/src/CaloShowerROMaker_module.cc +++ b/CaloMC/src/CaloShowerROMaker_module.cc @@ -212,7 +212,7 @@ namespace mu2e { float mbtime = GlobalConstantsHandle()->getNominalDRPeriod(); const Calorimeter& cal = *(GeomHandle()); - const float crystalLength = cal.caloInfo().getDouble("crystalZLength"); + const float crystalLength = cal.G4Info().get("crystalZLength"); std::map> simEntriesMap; diagSummary diagSum; diff --git a/CaloMC/src/CaloShowerStepMaker_module.cc b/CaloMC/src/CaloShowerStepMaker_module.cc index 394506e31f..be5513302e 100644 --- a/CaloMC/src/CaloShowerStepMaker_module.cc +++ b/CaloMC/src/CaloShowerStepMaker_module.cc @@ -119,7 +119,6 @@ namespace mu2e { void makeCompressedHits (const HandleVector&, CaloShowerStepCollection&, SimParticlePtrCollection&); void collectStepBySimAncestor (const Calorimeter&, const PhysicalVolumeMultiHelper&, const HandleVector&, std::map&); void collectStepBySim (const HandleVector&, SimStepMap&); - bool isInsideCalorimeter (const Calorimeter& cal, const PhysicalVolumeMultiHelper&, const SimPtr&); void compressSteps (const Calorimeter&, CaloShowerStepCollection&, int, const SimPtr&, std::vector&); void fillHisto1 (const Calorimeter&, const SimPtr&, const std::set&); void fillHisto2 (int, float, const SimPtr&); @@ -248,7 +247,7 @@ namespace mu2e { PhysicalVolumeMultiHelper vi(vols_); const Calorimeter& cal = *(GeomHandle()); - zSliceSize_ = cal.caloInfo().getDouble("crystalZLength")/float(numZSlices_)+1e-5; + zSliceSize_ = cal.G4Info().get("crystalZLength")/float(numZSlices_)+1e-5; //----------------------------------------------------------------- @@ -352,10 +351,10 @@ namespace mu2e { SimPtr sim = step.simParticle(); SimParticlePtrCollection inspectedSims; - while (sim->hasParent() && isInsideCalorimeter(cal, vi, sim) ) + while (sim->hasParent() && cal.isInsideAnyCrystal(sim->startPosition()) ) { //simparticle starting in one section and ending in another one see note above - if (!cal.geomUtil().isContainedSection(sim->startPosition(),sim->endPosition()) ) break; + if (!cal.isInsideSameDisk(sim->startPosition(),sim->endPosition()) ) break; const auto alreadyInspected = simToAncestorMap.find(sim); if (alreadyInspected != simToAncestorMap.end()) {sim = alreadyInspected->second; break;} @@ -374,17 +373,6 @@ namespace mu2e { } - - - - //------------------------------------------------------------------------------------------------------------------------- - bool CaloShowerStepMaker::isInsideCalorimeter(const Calorimeter& cal, const PhysicalVolumeMultiHelper& vi, - const art::Ptr& thisSimPtr) - { - if (usePhysVol_) return mapPhysVol_.find(&vi.startVolume(*thisSimPtr)) != mapPhysVol_.end(); - return cal.geomUtil().isInsideCalorimeter(thisSimPtr->startPosition()); - } - //----------------------------------------------------------------------------------------------------------------------------------------------- void CaloShowerStepMaker::collectStepBySim(const HandleVector& stepsHandles, std::map>& simStepMap) @@ -409,7 +397,7 @@ namespace mu2e { for (const StepPointMC* step : steps) { - CLHEP::Hep3Vector pos = cal.geomUtil().mu2eToCrystal(volId,step->position()); + CLHEP::Hep3Vector pos = cal.mu2eToCrystal(volId,step->position()); int idx = int(std::max(1e-6,pos.z())/zSliceSize_); if (buffer.entries(idx)>0 && (step->time()-buffer.t0(idx) > deltaTime_) ) @@ -443,8 +431,8 @@ namespace mu2e { //------------------------------------------------------------------------------------------------------------- void CaloShowerStepMaker::fillHisto1(const Calorimeter& cal, const art::Ptr& sim, const std::set>& infoSims) { - CLHEP::Hep3Vector startSection = cal.geomUtil().mu2eToDisk(0,sim->startPosition()); - CLHEP::Hep3Vector endSection = cal.geomUtil().mu2eToDisk(0,sim->endPosition()); + CLHEP::Hep3Vector startSection = cal.mu2eToDisk(0,sim->startPosition()); + CLHEP::Hep3Vector endSection = cal.mu2eToDisk(0,sim->endPosition()); double rStart = sqrt(startSection.x()*startSection.x()+startSection.y()*startSection.y()); double rEnd = sqrt(endSection.x()*endSection.x()+endSection.y()*endSection.y()); @@ -477,9 +465,9 @@ namespace mu2e { std::cout< imax_) throw cet::exception("Rethrow")<< "[CaloMC/ShowerStepUtil] Index out of bound " << i << std::endl; diff --git a/CaloReco/src/CaloHitMaker_module.cc b/CaloReco/src/CaloHitMaker_module.cc index 14af1d86b2..af86599008 100644 --- a/CaloReco/src/CaloHitMaker_module.cc +++ b/CaloReco/src/CaloHitMaker_module.cc @@ -113,7 +113,7 @@ namespace mu2e { // fill the map that associate for each crystal the corresponding CaloRecoDigi indexes - int nSiPM = cal.nCrystals()*cal.caloInfo().getInt("nSiPMPerCrystal"); + int nSiPM = cal.nCrystals()*cal.G4Info().get("nSiPMPerCrystal"); std::vector> hitMap(nSiPM,std::vector()); for (unsigned i=0; i< recoCaloDigis.size(); ++i) diff --git a/CalorimeterGeom/CMakeLists.txt b/CalorimeterGeom/CMakeLists.txt index 148cb2ad2b..25b8571ac2 100644 --- a/CalorimeterGeom/CMakeLists.txt +++ b/CalorimeterGeom/CMakeLists.txt @@ -1,13 +1,19 @@ cet_make_library( SOURCE - src/CaloGeomUtil.cc + src/CaloUtil.cc src/DiskCalorimeter.cc src/Disk.cc + src/DiskInfo.cc src/SquareMapper.cc src/SquareShiftMapper.cc LIBRARIES PUBLIC - + + Offline::ConfigTools + Offline::DataProducts Offline::Mu2eInterfaces + + CLHEP::CLHEP + cetlib_except::cetlib_except ) install_source(SUBDIRS src) diff --git a/CalorimeterGeom/data/crystalPos.txt b/CalorimeterGeom/data/crystalPos.txt new file mode 100644 index 0000000000..4263f129b4 --- /dev/null +++ b/CalorimeterGeom/data/crystalPos.txt @@ -0,0 +1,1348 @@ +0 0 -188.65 377.3 +0 1 -154.35 377.3 +0 2 -120.05 377.3 +0 3 120.05 377.3 +0 4 154.35 377.3 +0 5 188.65 377.3 +0 6 205.8 343 +0 7 205.8 -343 +0 8 188.65 -377.3 +0 9 154.35 -377.3 +0 10 120.05 -377.3 +0 11 -120.05 -377.3 +0 12 -154.35 -377.3 +0 13 -188.65 -377.3 +0 14 -205.8 -343 +0 15 -205.8 343 +0 16 -205.8 411.6 +0 17 -171.5 411.6 +0 18 -137.2 411.6 +0 19 -102.9 411.6 +0 20 -68.6 411.6 +0 21 -34.3 411.6 +0 22 0 411.6 +0 23 34.3 411.6 +0 24 68.6 411.6 +0 25 102.9 411.6 +0 26 137.2 411.6 +0 27 171.5 411.6 +0 28 205.8 411.6 +0 29 222.95 377.3 +0 30 240.1 343 +0 31 257.25 308.7 +0 32 394.45 34.3 +0 33 411.6 0 +0 34 394.45 -34.3 +0 35 257.25 -308.7 +0 36 240.1 -343 +0 37 222.95 -377.3 +0 38 205.8 -411.6 +0 39 171.5 -411.6 +0 40 137.2 -411.6 +0 41 102.9 -411.6 +0 42 68.6 -411.6 +0 43 34.3 -411.6 +0 44 0 -411.6 +0 45 -34.3 -411.6 +0 46 -68.6 -411.6 +0 47 -102.9 -411.6 +0 48 -137.2 -411.6 +0 49 -171.5 -411.6 +0 50 -205.8 -411.6 +0 51 -222.95 -377.3 +0 52 -240.1 -343 +0 53 -257.25 -308.7 +0 54 -394.45 -34.3 +0 55 -411.6 0 +0 56 -394.45 34.3 +0 57 -257.25 308.7 +0 58 -240.1 343 +0 59 -222.95 377.3 +0 60 -222.95 445.9 +0 61 -188.65 445.9 +0 62 -154.35 445.9 +0 63 -120.05 445.9 +0 64 -85.75 445.9 +0 65 -51.45 445.9 +0 66 -17.15 445.9 +0 67 17.15 445.9 +0 68 51.45 445.9 +0 69 85.75 445.9 +0 70 120.05 445.9 +0 71 154.35 445.9 +0 72 188.65 445.9 +0 73 222.95 445.9 +0 74 240.1 411.6 +0 75 257.25 377.3 +0 76 274.4 343 +0 77 291.55 308.7 +0 78 308.7 274.4 +0 79 325.85 240.1 +0 80 343 205.8 +0 81 360.15 171.5 +0 82 377.3 137.2 +0 83 394.45 102.9 +0 84 411.6 68.6 +0 85 428.75 34.3 +0 86 445.9 0 +0 87 428.75 -34.3 +0 88 411.6 -68.6 +0 89 394.45 -102.9 +0 90 377.3 -137.2 +0 91 360.15 -171.5 +0 92 343 -205.8 +0 93 325.85 -240.1 +0 94 308.7 -274.4 +0 95 291.55 -308.7 +0 96 274.4 -343 +0 97 257.25 -377.3 +0 98 240.1 -411.6 +0 99 222.95 -445.9 +0 100 188.65 -445.9 +0 101 154.35 -445.9 +0 102 120.05 -445.9 +0 103 85.75 -445.9 +0 104 51.45 -445.9 +0 105 17.15 -445.9 +0 106 -17.15 -445.9 +0 107 -51.45 -445.9 +0 108 -85.75 -445.9 +0 109 -120.05 -445.9 +0 110 -154.35 -445.9 +0 111 -188.65 -445.9 +0 112 -222.95 -445.9 +0 113 -240.1 -411.6 +0 114 -257.25 -377.3 +0 115 -274.4 -343 +0 116 -291.55 -308.7 +0 117 -308.7 -274.4 +0 118 -325.85 -240.1 +0 119 -343 -205.8 +0 120 -360.15 -171.5 +0 121 -377.3 -137.2 +0 122 -394.45 -102.9 +0 123 -411.6 -68.6 +0 124 -428.75 -34.3 +0 125 -445.9 0 +0 126 -428.75 34.3 +0 127 -411.6 68.6 +0 128 -394.45 102.9 +0 129 -377.3 137.2 +0 130 -360.15 171.5 +0 131 -343 205.8 +0 132 -325.85 240.1 +0 133 -308.7 274.4 +0 134 -291.55 308.7 +0 135 -274.4 343 +0 136 -257.25 377.3 +0 137 -240.1 411.6 +0 138 -240.1 480.2 +0 139 -205.8 480.2 +0 140 -171.5 480.2 +0 141 -137.2 480.2 +0 142 -102.9 480.2 +0 143 -68.6 480.2 +0 144 -34.3 480.2 +0 145 0 480.2 +0 146 34.3 480.2 +0 147 68.6 480.2 +0 148 102.9 480.2 +0 149 137.2 480.2 +0 150 171.5 480.2 +0 151 205.8 480.2 +0 152 240.1 480.2 +0 153 257.25 445.9 +0 154 274.4 411.6 +0 155 291.55 377.3 +0 156 308.7 343 +0 157 325.85 308.7 +0 158 343 274.4 +0 159 360.15 240.1 +0 160 377.3 205.8 +0 161 394.45 171.5 +0 162 411.6 137.2 +0 163 428.75 102.9 +0 164 445.9 68.6 +0 165 463.05 34.3 +0 166 480.2 0 +0 167 463.05 -34.3 +0 168 445.9 -68.6 +0 169 428.75 -102.9 +0 170 411.6 -137.2 +0 171 394.45 -171.5 +0 172 377.3 -205.8 +0 173 360.15 -240.1 +0 174 343 -274.4 +0 175 325.85 -308.7 +0 176 308.7 -343 +0 177 291.55 -377.3 +0 178 274.4 -411.6 +0 179 257.25 -445.9 +0 180 240.1 -480.2 +0 181 205.8 -480.2 +0 182 171.5 -480.2 +0 183 137.2 -480.2 +0 184 102.9 -480.2 +0 185 68.6 -480.2 +0 186 34.3 -480.2 +0 187 0 -480.2 +0 188 -34.3 -480.2 +0 189 -68.6 -480.2 +0 190 -102.9 -480.2 +0 191 -137.2 -480.2 +0 192 -171.5 -480.2 +0 193 -205.8 -480.2 +0 194 -240.1 -480.2 +0 195 -257.25 -445.9 +0 196 -274.4 -411.6 +0 197 -291.55 -377.3 +0 198 -308.7 -343 +0 199 -325.85 -308.7 +0 200 -343 -274.4 +0 201 -360.15 -240.1 +0 202 -377.3 -205.8 +0 203 -394.45 -171.5 +0 204 -411.6 -137.2 +0 205 -428.75 -102.9 +0 206 -445.9 -68.6 +0 207 -463.05 -34.3 +0 208 -480.2 0 +0 209 -463.05 34.3 +0 210 -445.9 68.6 +0 211 -428.75 102.9 +0 212 -411.6 137.2 +0 213 -394.45 171.5 +0 214 -377.3 205.8 +0 215 -360.15 240.1 +0 216 -343 274.4 +0 217 -325.85 308.7 +0 218 -308.7 343 +0 219 -291.55 377.3 +0 220 -274.4 411.6 +0 221 -257.25 445.9 +0 222 -257.25 514.5 +0 223 -222.95 514.5 +0 224 -188.65 514.5 +0 225 -154.35 514.5 +0 226 -120.05 514.5 +0 227 -85.75 514.5 +0 228 -51.45 514.5 +0 229 -17.15 514.5 +0 230 17.15 514.5 +0 231 51.45 514.5 +0 232 85.75 514.5 +0 233 120.05 514.5 +0 234 154.35 514.5 +0 235 188.65 514.5 +0 236 222.95 514.5 +0 237 257.25 514.5 +0 238 274.4 480.2 +0 239 291.55 445.9 +0 240 308.7 411.6 +0 241 325.85 377.3 +0 242 343 343 +0 243 360.15 308.7 +0 244 377.3 274.4 +0 245 394.45 240.1 +0 246 411.6 205.8 +0 247 428.75 171.5 +0 248 445.9 137.2 +0 249 463.05 102.9 +0 250 480.2 68.6 +0 251 497.35 34.3 +0 252 514.5 0 +0 253 497.35 -34.3 +0 254 480.2 -68.6 +0 255 463.05 -102.9 +0 256 445.9 -137.2 +0 257 428.75 -171.5 +0 258 411.6 -205.8 +0 259 394.45 -240.1 +0 260 377.3 -274.4 +0 261 360.15 -308.7 +0 262 343 -343 +0 263 325.85 -377.3 +0 264 308.7 -411.6 +0 265 291.55 -445.9 +0 266 274.4 -480.2 +0 267 257.25 -514.5 +0 268 222.95 -514.5 +0 269 188.65 -514.5 +0 270 154.35 -514.5 +0 271 120.05 -514.5 +0 272 85.75 -514.5 +0 273 51.45 -514.5 +0 274 17.15 -514.5 +0 275 -17.15 -514.5 +0 276 -51.45 -514.5 +0 277 -85.75 -514.5 +0 278 -120.05 -514.5 +0 279 -154.35 -514.5 +0 280 -188.65 -514.5 +0 281 -222.95 -514.5 +0 282 -257.25 -514.5 +0 283 -274.4 -480.2 +0 284 -291.55 -445.9 +0 285 -308.7 -411.6 +0 286 -325.85 -377.3 +0 287 -343 -343 +0 288 -360.15 -308.7 +0 289 -377.3 -274.4 +0 290 -394.45 -240.1 +0 291 -411.6 -205.8 +0 292 -428.75 -171.5 +0 293 -445.9 -137.2 +0 294 -463.05 -102.9 +0 295 -480.2 -68.6 +0 296 -497.35 -34.3 +0 297 -514.5 0 +0 298 -497.35 34.3 +0 299 -480.2 68.6 +0 300 -463.05 102.9 +0 301 -445.9 137.2 +0 302 -428.75 171.5 +0 303 -411.6 205.8 +0 304 -394.45 240.1 +0 305 -377.3 274.4 +0 306 -360.15 308.7 +0 307 -343 343 +0 308 -325.85 377.3 +0 309 -308.7 411.6 +0 310 -291.55 445.9 +0 311 -274.4 480.2 +0 312 -274.4 548.8 +0 313 -240.1 548.8 +0 314 -205.8 548.8 +0 315 -171.5 548.8 +0 316 -137.2 548.8 +0 317 -102.9 548.8 +0 318 -68.6 548.8 +0 319 -34.3 548.8 +0 320 0 548.8 +0 321 34.3 548.8 +0 322 68.6 548.8 +0 323 102.9 548.8 +0 324 137.2 548.8 +0 325 171.5 548.8 +0 326 205.8 548.8 +0 327 240.1 548.8 +0 328 274.4 548.8 +0 329 291.55 514.5 +0 330 308.7 480.2 +0 331 325.85 445.9 +0 332 343 411.6 +0 333 360.15 377.3 +0 334 377.3 343 +0 335 394.45 308.7 +0 336 411.6 274.4 +0 337 428.75 240.1 +0 338 445.9 205.8 +0 339 463.05 171.5 +0 340 480.2 137.2 +0 341 497.35 102.9 +0 342 514.5 68.6 +0 343 531.65 34.3 +0 344 548.8 0 +0 345 531.65 -34.3 +0 346 514.5 -68.6 +0 347 497.35 -102.9 +0 348 480.2 -137.2 +0 349 463.05 -171.5 +0 350 445.9 -205.8 +0 351 428.75 -240.1 +0 352 411.6 -274.4 +0 353 394.45 -308.7 +0 354 377.3 -343 +0 355 360.15 -377.3 +0 356 343 -411.6 +0 357 325.85 -445.9 +0 358 308.7 -480.2 +0 359 291.55 -514.5 +0 360 274.4 -548.8 +0 361 240.1 -548.8 +0 362 205.8 -548.8 +0 363 171.5 -548.8 +0 364 137.2 -548.8 +0 365 102.9 -548.8 +0 366 68.6 -548.8 +0 367 34.3 -548.8 +0 368 0 -548.8 +0 369 -34.3 -548.8 +0 370 -68.6 -548.8 +0 371 -102.9 -548.8 +0 372 -137.2 -548.8 +0 373 -171.5 -548.8 +0 374 -205.8 -548.8 +0 375 -240.1 -548.8 +0 376 -274.4 -548.8 +0 377 -291.55 -514.5 +0 378 -308.7 -480.2 +0 379 -325.85 -445.9 +0 380 -343 -411.6 +0 381 -360.15 -377.3 +0 382 -377.3 -343 +0 383 -394.45 -308.7 +0 384 -411.6 -274.4 +0 385 -428.75 -240.1 +0 386 -445.9 -205.8 +0 387 -463.05 -171.5 +0 388 -480.2 -137.2 +0 389 -497.35 -102.9 +0 390 -514.5 -68.6 +0 391 -531.65 -34.3 +0 392 -548.8 0 +0 393 -531.65 34.3 +0 394 -514.5 68.6 +0 395 -497.35 102.9 +0 396 -480.2 137.2 +0 397 -463.05 171.5 +0 398 -445.9 205.8 +0 399 -428.75 240.1 +0 400 -411.6 274.4 +0 401 -394.45 308.7 +0 402 -377.3 343 +0 403 -360.15 377.3 +0 404 -343 411.6 +0 405 -325.85 445.9 +0 406 -308.7 480.2 +0 407 -291.55 514.5 +0 408 -222.95 583.1 +0 409 -188.65 583.1 +0 410 -154.35 583.1 +0 411 -120.05 583.1 +0 412 -85.75 583.1 +0 413 -51.45 583.1 +0 414 -17.15 583.1 +0 415 17.15 583.1 +0 416 51.45 583.1 +0 417 85.75 583.1 +0 418 120.05 583.1 +0 419 154.35 583.1 +0 420 188.65 583.1 +0 421 222.95 583.1 +0 422 308.7 548.8 +0 423 325.85 514.5 +0 424 343 480.2 +0 425 360.15 445.9 +0 426 377.3 411.6 +0 427 394.45 377.3 +0 428 411.6 343 +0 429 428.75 308.7 +0 430 445.9 274.4 +0 431 463.05 240.1 +0 432 480.2 205.8 +0 433 497.35 171.5 +0 434 514.5 137.2 +0 435 531.65 102.9 +0 436 548.8 68.6 +0 437 565.95 34.3 +0 438 583.1 0 +0 439 565.95 -34.3 +0 440 548.8 -68.6 +0 441 531.65 -102.9 +0 442 514.5 -137.2 +0 443 497.35 -171.5 +0 444 480.2 -205.8 +0 445 463.05 -240.1 +0 446 445.9 -274.4 +0 447 428.75 -308.7 +0 448 411.6 -343 +0 449 394.45 -377.3 +0 450 377.3 -411.6 +0 451 360.15 -445.9 +0 452 343 -480.2 +0 453 325.85 -514.5 +0 454 308.7 -548.8 +0 455 222.95 -583.1 +0 456 188.65 -583.1 +0 457 154.35 -583.1 +0 458 120.05 -583.1 +0 459 85.75 -583.1 +0 460 51.45 -583.1 +0 461 17.15 -583.1 +0 462 -17.15 -583.1 +0 463 -51.45 -583.1 +0 464 -85.75 -583.1 +0 465 -120.05 -583.1 +0 466 -154.35 -583.1 +0 467 -188.65 -583.1 +0 468 -222.95 -583.1 +0 469 -308.7 -548.8 +0 470 -325.85 -514.5 +0 471 -343 -480.2 +0 472 -360.15 -445.9 +0 473 -377.3 -411.6 +0 474 -394.45 -377.3 +0 475 -411.6 -343 +0 476 -428.75 -308.7 +0 477 -445.9 -274.4 +0 478 -463.05 -240.1 +0 479 -480.2 -205.8 +0 480 -497.35 -171.5 +0 481 -514.5 -137.2 +0 482 -531.65 -102.9 +0 483 -548.8 -68.6 +0 484 -565.95 -34.3 +0 485 -583.1 0 +0 486 -565.95 34.3 +0 487 -548.8 68.6 +0 488 -531.65 102.9 +0 489 -514.5 137.2 +0 490 -497.35 171.5 +0 491 -480.2 205.8 +0 492 -463.05 240.1 +0 493 -445.9 274.4 +0 494 -428.75 308.7 +0 495 -411.6 343 +0 496 -394.45 377.3 +0 497 -377.3 411.6 +0 498 -360.15 445.9 +0 499 -343 480.2 +0 500 -325.85 514.5 +0 501 -308.7 548.8 +0 502 -137.2 617.4 +0 503 -102.9 617.4 +0 504 -68.6 617.4 +0 505 -34.3 617.4 +0 506 0 617.4 +0 507 34.3 617.4 +0 508 68.6 617.4 +0 509 102.9 617.4 +0 510 137.2 617.4 +0 511 360.15 514.5 +0 512 377.3 480.2 +0 513 394.45 445.9 +0 514 411.6 411.6 +0 515 428.75 377.3 +0 516 445.9 343 +0 517 463.05 308.7 +0 518 480.2 274.4 +0 519 497.35 240.1 +0 520 514.5 205.8 +0 521 531.65 171.5 +0 522 548.8 137.2 +0 523 565.95 102.9 +0 524 583.1 68.6 +0 525 600.25 34.3 +0 526 617.4 0 +0 527 600.25 -34.3 +0 528 583.1 -68.6 +0 529 565.95 -102.9 +0 530 548.8 -137.2 +0 531 531.65 -171.5 +0 532 514.5 -205.8 +0 533 497.35 -240.1 +0 534 480.2 -274.4 +0 535 463.05 -308.7 +0 536 445.9 -343 +0 537 428.75 -377.3 +0 538 411.6 -411.6 +0 539 394.45 -445.9 +0 540 377.3 -480.2 +0 541 360.15 -514.5 +0 542 137.2 -617.4 +0 543 102.9 -617.4 +0 544 68.6 -617.4 +0 545 34.3 -617.4 +0 546 0 -617.4 +0 547 -34.3 -617.4 +0 548 -68.6 -617.4 +0 549 -102.9 -617.4 +0 550 -137.2 -617.4 +0 551 -360.15 -514.5 +0 552 -377.3 -480.2 +0 553 -394.45 -445.9 +0 554 -411.6 -411.6 +0 555 -428.75 -377.3 +0 556 -445.9 -343 +0 557 -463.05 -308.7 +0 558 -480.2 -274.4 +0 559 -497.35 -240.1 +0 560 -514.5 -205.8 +0 561 -531.65 -171.5 +0 562 -548.8 -137.2 +0 563 -565.95 -102.9 +0 564 -583.1 -68.6 +0 565 -600.25 -34.3 +0 566 -617.4 0 +0 567 -600.25 34.3 +0 568 -583.1 68.6 +0 569 -565.95 102.9 +0 570 -548.8 137.2 +0 571 -531.65 171.5 +0 572 -514.5 205.8 +0 573 -497.35 240.1 +0 574 -480.2 274.4 +0 575 -463.05 308.7 +0 576 -445.9 343 +0 577 -428.75 377.3 +0 578 -411.6 411.6 +0 579 -394.45 445.9 +0 580 -377.3 480.2 +0 581 -360.15 514.5 +0 582 411.6 480.2 +0 583 428.75 445.9 +0 584 445.9 411.6 +0 585 463.05 377.3 +0 586 480.2 343 +0 587 497.35 308.7 +0 588 514.5 274.4 +0 589 531.65 240.1 +0 590 548.8 205.8 +0 591 565.95 171.5 +0 592 583.1 137.2 +0 593 600.25 102.9 +0 594 617.4 68.6 +0 595 634.55 34.3 +0 596 634.55 -34.3 +0 597 617.4 -68.6 +0 598 600.25 -102.9 +0 599 583.1 -137.2 +0 600 565.95 -171.5 +0 601 548.8 -205.8 +0 602 531.65 -240.1 +0 603 514.5 -274.4 +0 604 497.35 -308.7 +0 605 480.2 -343 +0 606 463.05 -377.3 +0 607 445.9 -411.6 +0 608 428.75 -445.9 +0 609 411.6 -480.2 +0 610 -411.6 -480.2 +0 611 -428.75 -445.9 +0 612 -445.9 -411.6 +0 613 -463.05 -377.3 +0 614 -480.2 -343 +0 615 -497.35 -308.7 +0 616 -514.5 -274.4 +0 617 -531.65 -240.1 +0 618 -548.8 -205.8 +0 619 -565.95 -171.5 +0 620 -583.1 -137.2 +0 621 -600.25 -102.9 +0 622 -617.4 -68.6 +0 623 -634.55 -34.3 +0 624 -634.55 34.3 +0 625 -617.4 68.6 +0 626 -600.25 102.9 +0 627 -583.1 137.2 +0 628 -565.95 171.5 +0 629 -548.8 205.8 +0 630 -531.65 240.1 +0 631 -514.5 274.4 +0 632 -497.35 308.7 +0 633 -480.2 343 +0 634 -463.05 377.3 +0 635 -445.9 411.6 +0 636 -428.75 445.9 +0 637 -411.6 480.2 +0 638 480.2 411.6 +0 639 497.35 377.3 +0 640 514.5 343 +0 641 531.65 308.7 +0 642 548.8 274.4 +0 643 565.95 240.1 +0 644 583.1 205.8 +0 645 600.25 171.5 +0 646 617.4 137.2 +0 647 617.4 -137.2 +0 648 600.25 -171.5 +0 649 583.1 -205.8 +0 650 565.95 -240.1 +0 651 548.8 -274.4 +0 652 531.65 -308.7 +0 653 514.5 -343 +0 654 497.35 -377.3 +0 655 480.2 -411.6 +0 656 -480.2 -411.6 +0 657 -497.35 -377.3 +0 658 -514.5 -343 +0 659 -531.65 -308.7 +0 660 -548.8 -274.4 +0 661 -565.95 -240.1 +0 662 -583.1 -205.8 +0 663 -600.25 -171.5 +0 664 -617.4 -137.2 +0 665 -617.4 137.2 +0 666 -600.25 171.5 +0 667 -583.1 205.8 +0 668 -565.95 240.1 +0 669 -548.8 274.4 +0 670 -531.65 308.7 +0 671 -514.5 343 +0 672 -497.35 377.3 +0 673 -480.2 411.6 +1 0 -188.65 377.3 +1 1 -154.35 377.3 +1 2 -120.05 377.3 +1 3 120.05 377.3 +1 4 154.35 377.3 +1 5 188.65 377.3 +1 6 205.8 343 +1 7 205.8 -343 +1 8 188.65 -377.3 +1 9 154.35 -377.3 +1 10 120.05 -377.3 +1 11 -120.05 -377.3 +1 12 -154.35 -377.3 +1 13 -188.65 -377.3 +1 14 -205.8 -343 +1 15 -205.8 343 +1 16 -205.8 411.6 +1 17 -171.5 411.6 +1 18 -137.2 411.6 +1 19 -102.9 411.6 +1 20 -68.6 411.6 +1 21 -34.3 411.6 +1 22 0 411.6 +1 23 34.3 411.6 +1 24 68.6 411.6 +1 25 102.9 411.6 +1 26 137.2 411.6 +1 27 171.5 411.6 +1 28 205.8 411.6 +1 29 222.95 377.3 +1 30 240.1 343 +1 31 257.25 308.7 +1 32 394.45 34.3 +1 33 411.6 0 +1 34 394.45 -34.3 +1 35 257.25 -308.7 +1 36 240.1 -343 +1 37 222.95 -377.3 +1 38 205.8 -411.6 +1 39 171.5 -411.6 +1 40 137.2 -411.6 +1 41 102.9 -411.6 +1 42 68.6 -411.6 +1 43 34.3 -411.6 +1 44 0 -411.6 +1 45 -34.3 -411.6 +1 46 -68.6 -411.6 +1 47 -102.9 -411.6 +1 48 -137.2 -411.6 +1 49 -171.5 -411.6 +1 50 -205.8 -411.6 +1 51 -222.95 -377.3 +1 52 -240.1 -343 +1 53 -257.25 -308.7 +1 54 -394.45 -34.3 +1 55 -411.6 0 +1 56 -394.45 34.3 +1 57 -257.25 308.7 +1 58 -240.1 343 +1 59 -222.95 377.3 +1 60 -222.95 445.9 +1 61 -188.65 445.9 +1 62 -154.35 445.9 +1 63 -120.05 445.9 +1 64 -85.75 445.9 +1 65 -51.45 445.9 +1 66 -17.15 445.9 +1 67 17.15 445.9 +1 68 51.45 445.9 +1 69 85.75 445.9 +1 70 120.05 445.9 +1 71 154.35 445.9 +1 72 188.65 445.9 +1 73 222.95 445.9 +1 74 240.1 411.6 +1 75 257.25 377.3 +1 76 274.4 343 +1 77 291.55 308.7 +1 78 308.7 274.4 +1 79 325.85 240.1 +1 80 343 205.8 +1 81 360.15 171.5 +1 82 377.3 137.2 +1 83 394.45 102.9 +1 84 411.6 68.6 +1 85 428.75 34.3 +1 86 445.9 0 +1 87 428.75 -34.3 +1 88 411.6 -68.6 +1 89 394.45 -102.9 +1 90 377.3 -137.2 +1 91 360.15 -171.5 +1 92 343 -205.8 +1 93 325.85 -240.1 +1 94 308.7 -274.4 +1 95 291.55 -308.7 +1 96 274.4 -343 +1 97 257.25 -377.3 +1 98 240.1 -411.6 +1 99 222.95 -445.9 +1 100 188.65 -445.9 +1 101 154.35 -445.9 +1 102 120.05 -445.9 +1 103 85.75 -445.9 +1 104 51.45 -445.9 +1 105 17.15 -445.9 +1 106 -17.15 -445.9 +1 107 -51.45 -445.9 +1 108 -85.75 -445.9 +1 109 -120.05 -445.9 +1 110 -154.35 -445.9 +1 111 -188.65 -445.9 +1 112 -222.95 -445.9 +1 113 -240.1 -411.6 +1 114 -257.25 -377.3 +1 115 -274.4 -343 +1 116 -291.55 -308.7 +1 117 -308.7 -274.4 +1 118 -325.85 -240.1 +1 119 -343 -205.8 +1 120 -360.15 -171.5 +1 121 -377.3 -137.2 +1 122 -394.45 -102.9 +1 123 -411.6 -68.6 +1 124 -428.75 -34.3 +1 125 -445.9 0 +1 126 -428.75 34.3 +1 127 -411.6 68.6 +1 128 -394.45 102.9 +1 129 -377.3 137.2 +1 130 -360.15 171.5 +1 131 -343 205.8 +1 132 -325.85 240.1 +1 133 -308.7 274.4 +1 134 -291.55 308.7 +1 135 -274.4 343 +1 136 -257.25 377.3 +1 137 -240.1 411.6 +1 138 -240.1 480.2 +1 139 -205.8 480.2 +1 140 -171.5 480.2 +1 141 -137.2 480.2 +1 142 -102.9 480.2 +1 143 -68.6 480.2 +1 144 -34.3 480.2 +1 145 0 480.2 +1 146 34.3 480.2 +1 147 68.6 480.2 +1 148 102.9 480.2 +1 149 137.2 480.2 +1 150 171.5 480.2 +1 151 205.8 480.2 +1 152 240.1 480.2 +1 153 257.25 445.9 +1 154 274.4 411.6 +1 155 291.55 377.3 +1 156 308.7 343 +1 157 325.85 308.7 +1 158 343 274.4 +1 159 360.15 240.1 +1 160 377.3 205.8 +1 161 394.45 171.5 +1 162 411.6 137.2 +1 163 428.75 102.9 +1 164 445.9 68.6 +1 165 463.05 34.3 +1 166 480.2 0 +1 167 463.05 -34.3 +1 168 445.9 -68.6 +1 169 428.75 -102.9 +1 170 411.6 -137.2 +1 171 394.45 -171.5 +1 172 377.3 -205.8 +1 173 360.15 -240.1 +1 174 343 -274.4 +1 175 325.85 -308.7 +1 176 308.7 -343 +1 177 291.55 -377.3 +1 178 274.4 -411.6 +1 179 257.25 -445.9 +1 180 240.1 -480.2 +1 181 205.8 -480.2 +1 182 171.5 -480.2 +1 183 137.2 -480.2 +1 184 102.9 -480.2 +1 185 68.6 -480.2 +1 186 34.3 -480.2 +1 187 0 -480.2 +1 188 -34.3 -480.2 +1 189 -68.6 -480.2 +1 190 -102.9 -480.2 +1 191 -137.2 -480.2 +1 192 -171.5 -480.2 +1 193 -205.8 -480.2 +1 194 -240.1 -480.2 +1 195 -257.25 -445.9 +1 196 -274.4 -411.6 +1 197 -291.55 -377.3 +1 198 -308.7 -343 +1 199 -325.85 -308.7 +1 200 -343 -274.4 +1 201 -360.15 -240.1 +1 202 -377.3 -205.8 +1 203 -394.45 -171.5 +1 204 -411.6 -137.2 +1 205 -428.75 -102.9 +1 206 -445.9 -68.6 +1 207 -463.05 -34.3 +1 208 -480.2 0 +1 209 -463.05 34.3 +1 210 -445.9 68.6 +1 211 -428.75 102.9 +1 212 -411.6 137.2 +1 213 -394.45 171.5 +1 214 -377.3 205.8 +1 215 -360.15 240.1 +1 216 -343 274.4 +1 217 -325.85 308.7 +1 218 -308.7 343 +1 219 -291.55 377.3 +1 220 -274.4 411.6 +1 221 -257.25 445.9 +1 222 -257.25 514.5 +1 223 -222.95 514.5 +1 224 -188.65 514.5 +1 225 -154.35 514.5 +1 226 -120.05 514.5 +1 227 -85.75 514.5 +1 228 -51.45 514.5 +1 229 -17.15 514.5 +1 230 17.15 514.5 +1 231 51.45 514.5 +1 232 85.75 514.5 +1 233 120.05 514.5 +1 234 154.35 514.5 +1 235 188.65 514.5 +1 236 222.95 514.5 +1 237 257.25 514.5 +1 238 274.4 480.2 +1 239 291.55 445.9 +1 240 308.7 411.6 +1 241 325.85 377.3 +1 242 343 343 +1 243 360.15 308.7 +1 244 377.3 274.4 +1 245 394.45 240.1 +1 246 411.6 205.8 +1 247 428.75 171.5 +1 248 445.9 137.2 +1 249 463.05 102.9 +1 250 480.2 68.6 +1 251 497.35 34.3 +1 252 514.5 0 +1 253 497.35 -34.3 +1 254 480.2 -68.6 +1 255 463.05 -102.9 +1 256 445.9 -137.2 +1 257 428.75 -171.5 +1 258 411.6 -205.8 +1 259 394.45 -240.1 +1 260 377.3 -274.4 +1 261 360.15 -308.7 +1 262 343 -343 +1 263 325.85 -377.3 +1 264 308.7 -411.6 +1 265 291.55 -445.9 +1 266 274.4 -480.2 +1 267 257.25 -514.5 +1 268 222.95 -514.5 +1 269 188.65 -514.5 +1 270 154.35 -514.5 +1 271 120.05 -514.5 +1 272 85.75 -514.5 +1 273 51.45 -514.5 +1 274 17.15 -514.5 +1 275 -17.15 -514.5 +1 276 -51.45 -514.5 +1 277 -85.75 -514.5 +1 278 -120.05 -514.5 +1 279 -154.35 -514.5 +1 280 -188.65 -514.5 +1 281 -222.95 -514.5 +1 282 -257.25 -514.5 +1 283 -274.4 -480.2 +1 284 -291.55 -445.9 +1 285 -308.7 -411.6 +1 286 -325.85 -377.3 +1 287 -343 -343 +1 288 -360.15 -308.7 +1 289 -377.3 -274.4 +1 290 -394.45 -240.1 +1 291 -411.6 -205.8 +1 292 -428.75 -171.5 +1 293 -445.9 -137.2 +1 294 -463.05 -102.9 +1 295 -480.2 -68.6 +1 296 -497.35 -34.3 +1 297 -514.5 0 +1 298 -497.35 34.3 +1 299 -480.2 68.6 +1 300 -463.05 102.9 +1 301 -445.9 137.2 +1 302 -428.75 171.5 +1 303 -411.6 205.8 +1 304 -394.45 240.1 +1 305 -377.3 274.4 +1 306 -360.15 308.7 +1 307 -343 343 +1 308 -325.85 377.3 +1 309 -308.7 411.6 +1 310 -291.55 445.9 +1 311 -274.4 480.2 +1 312 -274.4 548.8 +1 313 -240.1 548.8 +1 314 -205.8 548.8 +1 315 -171.5 548.8 +1 316 -137.2 548.8 +1 317 -102.9 548.8 +1 318 -68.6 548.8 +1 319 -34.3 548.8 +1 320 0 548.8 +1 321 34.3 548.8 +1 322 68.6 548.8 +1 323 102.9 548.8 +1 324 137.2 548.8 +1 325 171.5 548.8 +1 326 205.8 548.8 +1 327 240.1 548.8 +1 328 274.4 548.8 +1 329 291.55 514.5 +1 330 308.7 480.2 +1 331 325.85 445.9 +1 332 343 411.6 +1 333 360.15 377.3 +1 334 377.3 343 +1 335 394.45 308.7 +1 336 411.6 274.4 +1 337 428.75 240.1 +1 338 445.9 205.8 +1 339 463.05 171.5 +1 340 480.2 137.2 +1 341 497.35 102.9 +1 342 514.5 68.6 +1 343 531.65 34.3 +1 344 548.8 0 +1 345 531.65 -34.3 +1 346 514.5 -68.6 +1 347 497.35 -102.9 +1 348 480.2 -137.2 +1 349 463.05 -171.5 +1 350 445.9 -205.8 +1 351 428.75 -240.1 +1 352 411.6 -274.4 +1 353 394.45 -308.7 +1 354 377.3 -343 +1 355 360.15 -377.3 +1 356 343 -411.6 +1 357 325.85 -445.9 +1 358 308.7 -480.2 +1 359 291.55 -514.5 +1 360 274.4 -548.8 +1 361 240.1 -548.8 +1 362 205.8 -548.8 +1 363 171.5 -548.8 +1 364 137.2 -548.8 +1 365 102.9 -548.8 +1 366 68.6 -548.8 +1 367 34.3 -548.8 +1 368 0 -548.8 +1 369 -34.3 -548.8 +1 370 -68.6 -548.8 +1 371 -102.9 -548.8 +1 372 -137.2 -548.8 +1 373 -171.5 -548.8 +1 374 -205.8 -548.8 +1 375 -240.1 -548.8 +1 376 -274.4 -548.8 +1 377 -291.55 -514.5 +1 378 -308.7 -480.2 +1 379 -325.85 -445.9 +1 380 -343 -411.6 +1 381 -360.15 -377.3 +1 382 -377.3 -343 +1 383 -394.45 -308.7 +1 384 -411.6 -274.4 +1 385 -428.75 -240.1 +1 386 -445.9 -205.8 +1 387 -463.05 -171.5 +1 388 -480.2 -137.2 +1 389 -497.35 -102.9 +1 390 -514.5 -68.6 +1 391 -531.65 -34.3 +1 392 -548.8 0 +1 393 -531.65 34.3 +1 394 -514.5 68.6 +1 395 -497.35 102.9 +1 396 -480.2 137.2 +1 397 -463.05 171.5 +1 398 -445.9 205.8 +1 399 -428.75 240.1 +1 400 -411.6 274.4 +1 401 -394.45 308.7 +1 402 -377.3 343 +1 403 -360.15 377.3 +1 404 -343 411.6 +1 405 -325.85 445.9 +1 406 -308.7 480.2 +1 407 -291.55 514.5 +1 408 -222.95 583.1 +1 409 -188.65 583.1 +1 410 -154.35 583.1 +1 411 -120.05 583.1 +1 412 -85.75 583.1 +1 413 -51.45 583.1 +1 414 -17.15 583.1 +1 415 17.15 583.1 +1 416 51.45 583.1 +1 417 85.75 583.1 +1 418 120.05 583.1 +1 419 154.35 583.1 +1 420 188.65 583.1 +1 421 222.95 583.1 +1 422 308.7 548.8 +1 423 325.85 514.5 +1 424 343 480.2 +1 425 360.15 445.9 +1 426 377.3 411.6 +1 427 394.45 377.3 +1 428 411.6 343 +1 429 428.75 308.7 +1 430 445.9 274.4 +1 431 463.05 240.1 +1 432 480.2 205.8 +1 433 497.35 171.5 +1 434 514.5 137.2 +1 435 531.65 102.9 +1 436 548.8 68.6 +1 437 565.95 34.3 +1 438 583.1 0 +1 439 565.95 -34.3 +1 440 548.8 -68.6 +1 441 531.65 -102.9 +1 442 514.5 -137.2 +1 443 497.35 -171.5 +1 444 480.2 -205.8 +1 445 463.05 -240.1 +1 446 445.9 -274.4 +1 447 428.75 -308.7 +1 448 411.6 -343 +1 449 394.45 -377.3 +1 450 377.3 -411.6 +1 451 360.15 -445.9 +1 452 343 -480.2 +1 453 325.85 -514.5 +1 454 308.7 -548.8 +1 455 222.95 -583.1 +1 456 188.65 -583.1 +1 457 154.35 -583.1 +1 458 120.05 -583.1 +1 459 85.75 -583.1 +1 460 51.45 -583.1 +1 461 17.15 -583.1 +1 462 -17.15 -583.1 +1 463 -51.45 -583.1 +1 464 -85.75 -583.1 +1 465 -120.05 -583.1 +1 466 -154.35 -583.1 +1 467 -188.65 -583.1 +1 468 -222.95 -583.1 +1 469 -308.7 -548.8 +1 470 -325.85 -514.5 +1 471 -343 -480.2 +1 472 -360.15 -445.9 +1 473 -377.3 -411.6 +1 474 -394.45 -377.3 +1 475 -411.6 -343 +1 476 -428.75 -308.7 +1 477 -445.9 -274.4 +1 478 -463.05 -240.1 +1 479 -480.2 -205.8 +1 480 -497.35 -171.5 +1 481 -514.5 -137.2 +1 482 -531.65 -102.9 +1 483 -548.8 -68.6 +1 484 -565.95 -34.3 +1 485 -583.1 0 +1 486 -565.95 34.3 +1 487 -548.8 68.6 +1 488 -531.65 102.9 +1 489 -514.5 137.2 +1 490 -497.35 171.5 +1 491 -480.2 205.8 +1 492 -463.05 240.1 +1 493 -445.9 274.4 +1 494 -428.75 308.7 +1 495 -411.6 343 +1 496 -394.45 377.3 +1 497 -377.3 411.6 +1 498 -360.15 445.9 +1 499 -343 480.2 +1 500 -325.85 514.5 +1 501 -308.7 548.8 +1 502 -137.2 617.4 +1 503 -102.9 617.4 +1 504 -68.6 617.4 +1 505 -34.3 617.4 +1 506 0 617.4 +1 507 34.3 617.4 +1 508 68.6 617.4 +1 509 102.9 617.4 +1 510 137.2 617.4 +1 511 360.15 514.5 +1 512 377.3 480.2 +1 513 394.45 445.9 +1 514 411.6 411.6 +1 515 428.75 377.3 +1 516 445.9 343 +1 517 463.05 308.7 +1 518 480.2 274.4 +1 519 497.35 240.1 +1 520 514.5 205.8 +1 521 531.65 171.5 +1 522 548.8 137.2 +1 523 565.95 102.9 +1 524 583.1 68.6 +1 525 600.25 34.3 +1 526 617.4 0 +1 527 600.25 -34.3 +1 528 583.1 -68.6 +1 529 565.95 -102.9 +1 530 548.8 -137.2 +1 531 531.65 -171.5 +1 532 514.5 -205.8 +1 533 497.35 -240.1 +1 534 480.2 -274.4 +1 535 463.05 -308.7 +1 536 445.9 -343 +1 537 428.75 -377.3 +1 538 411.6 -411.6 +1 539 394.45 -445.9 +1 540 377.3 -480.2 +1 541 360.15 -514.5 +1 542 137.2 -617.4 +1 543 102.9 -617.4 +1 544 68.6 -617.4 +1 545 34.3 -617.4 +1 546 0 -617.4 +1 547 -34.3 -617.4 +1 548 -68.6 -617.4 +1 549 -102.9 -617.4 +1 550 -137.2 -617.4 +1 551 -360.15 -514.5 +1 552 -377.3 -480.2 +1 553 -394.45 -445.9 +1 554 -411.6 -411.6 +1 555 -428.75 -377.3 +1 556 -445.9 -343 +1 557 -463.05 -308.7 +1 558 -480.2 -274.4 +1 559 -497.35 -240.1 +1 560 -514.5 -205.8 +1 561 -531.65 -171.5 +1 562 -548.8 -137.2 +1 563 -565.95 -102.9 +1 564 -583.1 -68.6 +1 565 -600.25 -34.3 +1 566 -617.4 0 +1 567 -600.25 34.3 +1 568 -583.1 68.6 +1 569 -565.95 102.9 +1 570 -548.8 137.2 +1 571 -531.65 171.5 +1 572 -514.5 205.8 +1 573 -497.35 240.1 +1 574 -480.2 274.4 +1 575 -463.05 308.7 +1 576 -445.9 343 +1 577 -428.75 377.3 +1 578 -411.6 411.6 +1 579 -394.45 445.9 +1 580 -377.3 480.2 +1 581 -360.15 514.5 +1 582 411.6 480.2 +1 583 428.75 445.9 +1 584 445.9 411.6 +1 585 463.05 377.3 +1 586 480.2 343 +1 587 497.35 308.7 +1 588 514.5 274.4 +1 589 531.65 240.1 +1 590 548.8 205.8 +1 591 565.95 171.5 +1 592 583.1 137.2 +1 593 600.25 102.9 +1 594 617.4 68.6 +1 595 634.55 34.3 +1 596 634.55 -34.3 +1 597 617.4 -68.6 +1 598 600.25 -102.9 +1 599 583.1 -137.2 +1 600 565.95 -171.5 +1 601 548.8 -205.8 +1 602 531.65 -240.1 +1 603 514.5 -274.4 +1 604 497.35 -308.7 +1 605 480.2 -343 +1 606 463.05 -377.3 +1 607 445.9 -411.6 +1 608 428.75 -445.9 +1 609 411.6 -480.2 +1 610 -411.6 -480.2 +1 611 -428.75 -445.9 +1 612 -445.9 -411.6 +1 613 -463.05 -377.3 +1 614 -480.2 -343 +1 615 -497.35 -308.7 +1 616 -514.5 -274.4 +1 617 -531.65 -240.1 +1 618 -548.8 -205.8 +1 619 -565.95 -171.5 +1 620 -583.1 -137.2 +1 621 -600.25 -102.9 +1 622 -617.4 -68.6 +1 623 -634.55 -34.3 +1 624 -634.55 34.3 +1 625 -617.4 68.6 +1 626 -600.25 102.9 +1 627 -583.1 137.2 +1 628 -565.95 171.5 +1 629 -548.8 205.8 +1 630 -531.65 240.1 +1 631 -514.5 274.4 +1 632 -497.35 308.7 +1 633 -480.2 343 +1 634 -463.05 377.3 +1 635 -445.9 411.6 +1 636 -428.75 445.9 +1 637 -411.6 480.2 +1 638 480.2 411.6 +1 639 497.35 377.3 +1 640 514.5 343 +1 641 531.65 308.7 +1 642 548.8 274.4 +1 643 565.95 240.1 +1 644 583.1 205.8 +1 645 600.25 171.5 +1 646 617.4 137.2 +1 647 617.4 -137.2 +1 648 600.25 -171.5 +1 649 583.1 -205.8 +1 650 565.95 -240.1 +1 651 548.8 -274.4 +1 652 531.65 -308.7 +1 653 514.5 -343 +1 654 497.35 -377.3 +1 655 480.2 -411.6 +1 656 -480.2 -411.6 +1 657 -497.35 -377.3 +1 658 -514.5 -343 +1 659 -531.65 -308.7 +1 660 -548.8 -274.4 +1 661 -565.95 -240.1 +1 662 -583.1 -205.8 +1 663 -600.25 -171.5 +1 664 -617.4 -137.2 +1 665 -617.4 137.2 +1 666 -600.25 171.5 +1 667 -583.1 205.8 +1 668 -565.95 240.1 +1 669 -548.8 274.4 +1 670 -531.65 308.7 +1 671 -514.5 343 +1 672 -497.35 377.3 +1 673 -480.2 411.6 diff --git a/CalorimeterGeom/inc/CaloG4Info.hh b/CalorimeterGeom/inc/CaloG4Info.hh new file mode 100644 index 0000000000..8151c19975 --- /dev/null +++ b/CalorimeterGeom/inc/CaloG4Info.hh @@ -0,0 +1,65 @@ +// +// Contains data for the calorimeter class. The non-critical fields are stored in a single keyed +// store; a few performance-critical fields accessed throughout the code are cached separately. +// +// Access is typed: set/get with T one of the supported value types (see Value below). +// See wiki for list of valid keys. +// +// Original author B. Echenard +// +#ifndef CalorimeterGeom_CaloG4Info_hh +#define CalorimeterGeom_CaloG4Info_hh + +#include "cetlib_except/exception.h" +#include "CLHEP/Vector/ThreeVector.h" + +#include +#include +#include +#include +#include +#include + + +namespace mu2e { + + class CaloG4Info { + public: + CaloG4Info() = default; + + template + void set(const std::string& key, T value) + { + static_assert(isSupported, "CaloG4Info::set - unsupported value type"); + data_[key] = std::move(value); + } + + template + const T& get(const std::string& key) const + { + static_assert(isSupported, "CaloG4Info::get - unsupported value type"); + const auto iter = data_.find(key); + if (iter == data_.end()) + throw cet::exception("CaloG4Info") << " unknown element " << key << "\n"; + if (!std::holds_alternative(iter->second)) + throw cet::exception("CaloG4Info") << " element " << key << " requested with the wrong type\n"; + return std::get(iter->second); + } + + + private: + using Value = std::variant, std::vector, + std::string, CLHEP::Hep3Vector>; + + template struct IsAlt : std::false_type {}; + template + struct IsAlt> : std::disjunction...> {}; + + template static constexpr bool isSupported = IsAlt::value; + + std::map data_; + }; +} + +#endif diff --git a/CalorimeterGeom/inc/CaloGeomUtil.hh b/CalorimeterGeom/inc/CaloGeomUtil.hh deleted file mode 100644 index 649328d9e7..0000000000 --- a/CalorimeterGeom/inc/CaloGeomUtil.hh +++ /dev/null @@ -1,63 +0,0 @@ -// -// Contains geometry utilities: coordinates transformations and checks if you are inside disks -// FF is the abrevation of front face -// -// Original author B. Echenard -// - -#ifndef CalorimeterGeom_CaloGeomUtil_hh -#define CalorimeterGeom_CaloGeomUtil_hh - -#include "Offline/CalorimeterGeom/inc/CaloGeomUtil.hh" -#include "Offline/CalorimeterGeom/inc/Disk.hh" -#include "CLHEP/Vector/ThreeVector.h" -#include "cetlib_except/exception.h" - -#include -#include - - -namespace mu2e { - - - class CaloGeomUtil { - - public: - CaloGeomUtil(const std::vector>& disks, - const std::vector& fullCrystalList); - - void origin (const CLHEP::Hep3Vector& vec) {origin_ = vec;} - void trackerCenter (const CLHEP::Hep3Vector& vec) {trackerCenter_ = vec;} - void crystalZLength(double value) {crystalZLength_ = value;} - - const CLHEP::Hep3Vector& origin() const {return origin_;} - const CLHEP::Hep3Vector& trackerCenter() const {return trackerCenter_;} - const Disk& disk(int i) const {return *disks_.at(i);} - const double crystalZLength() const {return crystalZLength_;} - - - CLHEP::Hep3Vector mu2eToCrystal(int crystalId, const CLHEP::Hep3Vector& pos) const; - CLHEP::Hep3Vector mu2eToDisk (int diskId, const CLHEP::Hep3Vector& pos) const; - CLHEP::Hep3Vector mu2eToDiskFF (int diskId, const CLHEP::Hep3Vector& pos) const; - CLHEP::Hep3Vector mu2eToTracker(const CLHEP::Hep3Vector& pos) const; - - CLHEP::Hep3Vector crystalToMu2e(int crystalId, const CLHEP::Hep3Vector& pos) const; - CLHEP::Hep3Vector diskToMu2e (int diskId, const CLHEP::Hep3Vector& pos) const; - CLHEP::Hep3Vector diskFFToMu2e (int diskId, const CLHEP::Hep3Vector& pos) const; - CLHEP::Hep3Vector trackerToMu2e(const CLHEP::Hep3Vector& pos) const; - - bool isInsideCalorimeter (const CLHEP::Hep3Vector& pos) const; - bool isInsideSection (int iDisk, const CLHEP::Hep3Vector& pos) const; - bool isContainedSection (const CLHEP::Hep3Vector& front, const CLHEP::Hep3Vector& back) const; - - - private: - const std::vector>& disks_; - const std::vector& fullCrystalList_; //non-owning pointers - CLHEP::Hep3Vector origin_; - CLHEP::Hep3Vector trackerCenter_; - double crystalZLength_; - }; -} - -#endif diff --git a/CalorimeterGeom/inc/CaloInfo.hh b/CalorimeterGeom/inc/CaloInfo.hh deleted file mode 100644 index e4f95b5173..0000000000 --- a/CalorimeterGeom/inc/CaloInfo.hh +++ /dev/null @@ -1,66 +0,0 @@ -// -// Contains data for the calorimeter class. The non-critical fields are saved into maps, and a few performance critical fields -// accessed throughout the code are cached for efficiency. -// -// See wiki for list of valid keys -// -// Original author B. Echenard -// -#ifndef CalorimeterGeom_CaloInfo_hh -#define CalorimeterGeom_CaloInfo_hh - -#include "cetlib_except/exception.h" -#include -#include -#include - - -namespace mu2e { - - template class CaloInfoData - { - public: - CaloInfoData() : data_() {}; - - void set(const std::string& key, const T& value) {data_[key] = value;} - - const T& get(const std::string& key) const - { - auto iter = data_.find(key); - if (iter == data_.end()) throw cet::exception("CaloInfo") << " unknown element "<second; - }; - - private: - std::map data_; - }; - - - - class CaloInfo { - - public: - CaloInfo() : dataBool_(), dataInt_(), dataDouble_(), dataVInt_(), dataVDouble_() {} - - void set(const std::string& key, bool value) {dataBool_.set(key,value);} - void set(const std::string& key, int value) {dataInt_.set(key,value);} - void set(const std::string& key, double value) {dataDouble_.set(key,value);} - void set(const std::string& key, const std::vector& value) {dataVInt_.set(key,value);} - void set(const std::string& key, const std::vector& value) {dataVDouble_.set(key,value);} - - const bool getBool (const std::string& key) const {return dataBool_.get(key);} - const int getInt (const std::string& key) const {return dataInt_.get(key);} - const double getDouble (const std::string& key) const {return dataDouble_.get(key);} - const std::vector getVInt (const std::string& key) const {return dataVInt_.get(key);} - const std::vector getVDouble(const std::string& key) const {return dataVDouble_.get(key);} - - private: - CaloInfoData dataBool_; - CaloInfoData dataInt_; - CaloInfoData dataDouble_; - CaloInfoData> dataVInt_; - CaloInfoData> dataVDouble_; - }; -} - -#endif diff --git a/CalorimeterGeom/inc/Calorimeter.hh b/CalorimeterGeom/inc/Calorimeter.hh index c00627f273..b4df32f937 100644 --- a/CalorimeterGeom/inc/Calorimeter.hh +++ b/CalorimeterGeom/inc/Calorimeter.hh @@ -1,54 +1,60 @@ +#ifndef CalorimeterGeom_Calorimeter_hh +#define CalorimeterGeom_Calorimeter_hh // // Interface to the disk calorimeter // Original author B. Echenard // - -#ifndef CalorimeterGeom_Calorimeter_hh -#define CalorimeterGeom_Calorimeter_hh - #include "Offline/Mu2eInterfaces/inc/Detector.hh" -#include "Offline/CalorimeterGeom/inc/Calorimeter.hh" -#include "Offline/CalorimeterGeom/inc/CaloGeomUtil.hh" -#include "Offline/CalorimeterGeom/inc/CaloInfo.hh" +#include "Offline/CalorimeterGeom/inc/CaloG4Info.hh" #include "Offline/CalorimeterGeom/inc/Disk.hh" #include "Offline/CalorimeterGeom/inc/Crystal.hh" #include "CLHEP/Vector/ThreeVector.h" #include +#include namespace mu2e { - class Calorimeter: virtual public Detector { + class Calorimeter: virtual public Detector { - using DiskPtr = std::shared_ptr; - using DiskPtrs = std::vector; - using Crystals = std::vector; - using CrystalPtrs = std::vector; + public: + using Disks = std::vector; + using Crystals = std::vector; + using Hep3Vector = CLHEP::Hep3Vector; - public: + constexpr static const char* cxname = {"Calorimeter"}; - //no constructor for this interface - virtual ~Calorimeter() = default; + //no constructor for this interface + virtual ~Calorimeter() = default; - virtual size_t nDisks() const = 0; - virtual const Disk& disk(size_t i) const = 0; - virtual const DiskPtrs& diskPtrs() const = 0; + virtual size_t nDisks() const = 0; + virtual const Disk& disk(unsigned i) const = 0; + virtual const Disks& disks() const = 0; - virtual size_t nCrystals() const = 0; - virtual const Crystal& crystal(size_t i) const = 0; - virtual const CrystalPtrs& crystalPtrs() const = 0; + virtual size_t nCrystals() const = 0; + virtual const Crystal& crystal(unsigned i) const = 0; + virtual const Crystals& crystals() const = 0; - virtual const CaloInfo& caloInfo() const = 0; - virtual const CaloGeomUtil& geomUtil() const = 0; + virtual const CaloG4Info& G4Info() const = 0; + virtual void print(std::ostream &os = std::cout) const = 0; - virtual const std::vector& neighbors(int crystalId) const = 0; - virtual const std::vector& nextNeighbors(int crystalId) const = 0; - virtual std::vector neighborsByLevel(int crystalId, int level) const = 0; - virtual int crystalIdxFromPosition(const CLHEP::Hep3Vector& pos) const = 0; - virtual int nearestIdxFromPosition(const CLHEP::Hep3Vector& pos) const = 0; + virtual std::vector neighbors(unsigned globalId, unsigned level) const = 0; + virtual bool isInsideAnyCrystal(const Hep3Vector& pos) const = 0; + virtual bool isInsideAnyDisk (const Hep3Vector& pos) const = 0; + virtual bool isInsideSameDisk (const Hep3Vector& front, + const Hep3Vector& back) const = 0; - virtual void print(std::ostream &os = std::cout) const = 0; - }; -} + virtual Hep3Vector mu2eToDisk (unsigned diskId, const Hep3Vector& p) const = 0; + virtual Hep3Vector diskToMu2e (unsigned diskId, const Hep3Vector& p) const = 0; + virtual Hep3Vector mu2eToDiskFF (unsigned diskId, const Hep3Vector& p) const = 0; + virtual Hep3Vector diskFFToMu2e (unsigned diskId, const Hep3Vector& p) const = 0; + virtual Hep3Vector mu2eToCrystal(unsigned crystalId, const Hep3Vector& p) const = 0; + virtual Hep3Vector crystalToMu2e(unsigned crystalId, const Hep3Vector& p) const = 0; + // These should REALLY go to the tracker + virtual Hep3Vector mu2eToTracker(const Hep3Vector& p) const = 0; + virtual Hep3Vector trackerToMu2e(const Hep3Vector& p) const = 0; + + }; +} #endif diff --git a/CalorimeterGeom/inc/Crystal.hh b/CalorimeterGeom/inc/Crystal.hh index 3101622d25..35a6dc1da2 100644 --- a/CalorimeterGeom/inc/Crystal.hh +++ b/CalorimeterGeom/inc/Crystal.hh @@ -1,68 +1,61 @@ +#ifndef CalorimeterGeom_Crystal_hh +#define CalorimeterGeom_Crystal_hh // // Hold information about a crystal -// Neighbors and position are given in the "Mu2e" frame, i.e the global frame -// localID and localPosition are given in the disk frame (i.e. local frame) -// Front face is sometimes abreviated FF (yes, I am lazy) +// ID, Neighbors, and position are given in the "Mu2e" global reference frame +// localID and localPosition are given in the local disk refernce frame // // Original author B Echenard // - -#ifndef CalorimeterGeom_Crystal_hh -#define CalorimeterGeom_Crystal_hh - +#include "Offline/DataProducts/inc/CaloConst.hh" #include "CLHEP/Vector/ThreeVector.h" #include namespace mu2e { - class Crystal { - - public: - Crystal(int localID, int diskID, const CLHEP::Hep3Vector& localPosition, - const CLHEP::Hep3Vector& idealLocalPosition, const CLHEP::Hep3Vector& size) : - localID_ (localID), - diskID_ (diskID), - size_ (size), - localPosition_ (localPosition), - idealLocalPosition_ (idealLocalPosition), - position_ (), - neighborsGlobal_ (), - nextNeighborsGlobal_ (), - neighborsGlobalRaw_ (), - nextNeighborsGlobalRaw_() - {} - - - int localID () const {return localID_;} - int diskID () const {return diskID_;} - const CLHEP::Hep3Vector& size () const {return size_;} - const CLHEP::Hep3Vector& localPosition() const {return localPosition_;} - const CLHEP::Hep3Vector& idealLocalPosition() const {return localPosition_;} - const CLHEP::Hep3Vector& position () const {return position_;} - const std::vector& neighbors () const {return neighborsGlobal_;} - const std::vector& nextNeighbors() const {return nextNeighborsGlobal_;} - - void setSize (const CLHEP::Hep3Vector& size) {size_ = size;} - void setLocalPosition(const CLHEP::Hep3Vector& pos) {localPosition_ = pos;} - void setPosition (const CLHEP::Hep3Vector& pos) {position_ = pos;} - void setNeighbors (const std::vector& list) {neighborsGlobal_=list;} - void setNextNeighbors(const std::vector& list) {nextNeighborsGlobal_=list;} - - - private: - int localID_; - int diskID_; - CLHEP::Hep3Vector size_; - CLHEP::Hep3Vector localPosition_; - CLHEP::Hep3Vector idealLocalPosition_; - CLHEP::Hep3Vector position_; - std::vector neighborsGlobal_; - std::vector nextNeighborsGlobal_; - std::vector neighborsGlobalRaw_; - std::vector nextNeighborsGlobalRaw_; - }; - + class Crystal { + public: + static constexpr unsigned invalidID_ = CaloConst::_invalid; + + Crystal(unsigned localID, unsigned diskID, + const CLHEP::Hep3Vector& localPosition, + const CLHEP::Hep3Vector& size) : + ID_(invalidID_), + localID_(localID), + diskID_(diskID), + size_(size), + localPosition_(localPosition), + position_(), + neighbors_(), + nextNeighbors_() + {} + + unsigned ID () const {return ID_;} + unsigned localID () const {return localID_;} + unsigned diskID () const {return diskID_;} + const CLHEP::Hep3Vector& size () const {return size_;} + const CLHEP::Hep3Vector& localPosition() const {return localPosition_;} + const CLHEP::Hep3Vector& position () const {return position_;} + const std::vector& neighbors () const {return neighbors_;} + const std::vector& nextNeighbors() const {return nextNeighbors_;} + + void setID (unsigned ID) {ID_ = ID;} + void setLocalID (unsigned localID) {localID_ = localID;} + void setDiskID (unsigned diskID) {diskID_ = diskID;} + void setLocalPosition(const CLHEP::Hep3Vector& pos) {localPosition_ = pos;} + void setPosition (const CLHEP::Hep3Vector& pos) {position_ = pos;} + void setNeighbors (const std::vector& list) {neighbors_ = list;} + void setNextNeighbors(const std::vector& list) {nextNeighbors_ = list;} + + private: + unsigned ID_; + unsigned localID_; + unsigned diskID_; + CLHEP::Hep3Vector size_; + CLHEP::Hep3Vector localPosition_; + CLHEP::Hep3Vector position_; + std::vector neighbors_; + std::vector nextNeighbors_; + }; } - - #endif diff --git a/CalorimeterGeom/inc/CrystalMapper.hh b/CalorimeterGeom/inc/CrystalMapper.hh index 1104a908df..8674f81389 100644 --- a/CalorimeterGeom/inc/CrystalMapper.hh +++ b/CalorimeterGeom/inc/CrystalMapper.hh @@ -3,8 +3,7 @@ // // Interface for classes describing the layout of the crystals in the disk. // -// Crystals could be square, hexagonal or triangular. -// Need to describe tesselation and shape +// Crystals could be square, hexagonal or triangular // // Original author B. Echenard // @@ -15,30 +14,20 @@ namespace mu2e { - class CrystalMapper { - - public: - //no constructor for this interface - virtual ~CrystalMapper() = default; - - virtual int nCrystalMax (int maxRing) const = 0; - - virtual CLHEP::Hep2Vector xyFromIndex (int thisIndex) const = 0; - virtual int indexFromXY (double x, double y) const = 0; - virtual int indexFromRowCol(int nRow, int nCol) const = 0; - virtual int rowFromIndex (int thisIndex) const = 0; - virtual int colFromIndex (int thisIndex) const = 0; - - virtual int numNeighbors (int level) const = 0; - virtual std::vector neighbors (int thisIndex, int level) const = 0; - - virtual const std::vector& apexX () const = 0; - virtual const std::vector& apexY () const = 0; - - virtual bool isInsideCrystal(double x, double y, - const CLHEP::Hep3Vector& pos, - const CLHEP::Hep3Vector& size) const = 0; - }; + class CrystalMapper { + public: + //no constructor for this interface + virtual ~CrystalMapper() = default; + + virtual unsigned nCrystalMax (unsigned maxRing) const = 0; + virtual CLHEP::Hep2Vector xyFromIndex (unsigned thisIndex) const = 0; + virtual unsigned indexFromXY (double x, double y) const = 0; + virtual unsigned indexFromRowCol(int nRow, int nCol) const = 0; + virtual int rowFromIndex (unsigned thisIndex) const = 0; + virtual int colFromIndex (unsigned thisIndex) const = 0; + virtual unsigned numNeighbors (unsigned level) const = 0; + virtual std::vector neighbors (unsigned thisIndex, unsigned level) const = 0; + }; } #endif diff --git a/CalorimeterGeom/inc/Disk.hh b/CalorimeterGeom/inc/Disk.hh index e92894e6b9..a8b3308567 100644 --- a/CalorimeterGeom/inc/Disk.hh +++ b/CalorimeterGeom/inc/Disk.hh @@ -1,12 +1,14 @@ -#ifndef CalorimeterGeom_CaloSection_hh -#define CalorimeterGeom_CaloSection_hh +#ifndef CalorimeterGeom_Disk_hh +#define CalorimeterGeom_Disk_hh // // Hold information about a disk in the calorimter. // +// The crystal numbering scheme start at the center of the disk +// to facilitate navigation +// // Original author B Echenard // - -#include "Offline/CalorimeterGeom/inc/DiskGeomInfo.hh" +#include "Offline/CalorimeterGeom/inc/DiskInfo.hh" #include "Offline/CalorimeterGeom/inc/Crystal.hh" #include "Offline/CalorimeterGeom/inc/CrystalMapper.hh" @@ -15,58 +17,66 @@ #include "CLHEP/Vector/Rotation.h" #include #include +#include +#include namespace mu2e { - class Disk { + class Disk { + public: + Disk(unsigned id, double rCrystalIn, double rCrystalOut, + double nominalCellSize, double nominalCellLength, unsigned offset, + const CLHEP::Hep3Vector& diskOriginToCrystalOrigin, + const std::string& crystalFileName); - public: - Disk(int id, double rCrystalIn, double rCrystalOut, double nominalCellSize, - double nominalCellLength, int offset, const CLHEP::Hep3Vector& diskOriginToCrystalOrigin); + unsigned id() const {return id_;} - int id() const {return id_;} + size_t nCrystals() const {return crystals_.size();} + const Crystal& crystal(unsigned i) const {return crystals_.at(i);} + Crystal& crystal(unsigned i) {return crystals_.at(i);} + unsigned crystalOffset() const {return crystalOffset_;} - size_t nCrystals() const {return crystalList_.size();} - int crystalOffset() const {return globalCrystalOffset_;} - const Crystal& crystal(int i) const {return crystalList_.at(i);} - Crystal& crystal(int i) {return crystalList_.at(i);} + const DiskInfo& diskInfo() const {return diskInfo_;} + DiskInfo& diskInfo() {return diskInfo_;} - const DiskGeomInfo& geomInfo() const {return geomInfo_;} - DiskGeomInfo& geomInfo() {return geomInfo_;} + std::vector neighbors (unsigned localId, unsigned level=1) const; + std::vector idxFromRow (int thisRow) const; + unsigned idxFromPosition(const CLHEP::Hep3Vector& pos) const; + bool isInsideDisk (const CLHEP::Hep3Vector& pos) const; + bool isInsideCrystal(const CLHEP::Hep3Vector& pos) const; - int idxFromPosition (double x, double y) const; - std::vector nearestNeighborsFromPos(double x, double y) const; - std::vector findLocalNeighbors (int crystalId, int level) const; - int idMinCrystalInside (int row) const; - int idMaxCrystalInside (int row) const; - void boundingBoxes (int thisRow,std::vector& params) const; + CLHEP::Hep3Vector mu2eToDisk (const CLHEP::Hep3Vector& p) const; + CLHEP::Hep3Vector diskToMu2e (const CLHEP::Hep3Vector& p) const; + CLHEP::Hep3Vector mu2eToDiskFF (const CLHEP::Hep3Vector& p) const; + CLHEP::Hep3Vector diskFFToMu2e (const CLHEP::Hep3Vector& p) const; + CLHEP::Hep3Vector mu2eToCrystal (unsigned localId, const CLHEP::Hep3Vector& p) const; + CLHEP::Hep3Vector crystalToMu2e (unsigned localId, const CLHEP::Hep3Vector& p) const; - void print(std::ostream& os = std::cout) const; + void moveCrystal (unsigned localId, const CLHEP::Hep3Vector& disp); + void moveDisk (const CLHEP::Hep3Vector& disp, + const CLHEP::HepRotation& rotation); - private: - void fillCrystalsIdeal (const CLHEP::Hep3Vector& crystalOriginInDisk, - double nominalCellLength); - void fillCrystals (const CLHEP::Hep3Vector& crystalOriginInDisk); - void fixCrystalPosition(); - bool isInsideDisk (double x, double y, double widthX, double widthY) const; - bool isInsideCrystal (int icry, double x, double y) const; - const bool isCrystalIdxValid (int i) const; - const bool isMapIdxValid (int i) const; + void print (std::ostream& os = std::cout) const; - int id_; - std::vector crystalList_; - DiskGeomInfo geomInfo_; - double radiusInCrystal_; - double radiusOutCrystal_; - double nominalCellSize_; - int globalCrystalOffset_; - std::shared_ptr crystalMap_; - std::vector mapToCrystal_; - std::vector crystalToMap_; - int rowMax_; - static constexpr int invalidID_ = -1; - }; -} + private: + void fillCrystals (const CLHEP::Hep3Vector&, double nominalCellLength, + const std::string& filename); + void checkPosition (const CLHEP::Hep3Vector& pos, const CLHEP::Hep3Vector& size, + unsigned crystalID) const; + bool isInsideCrystal (unsigned icry, const CLHEP::Hep3Vector& pos) const; + bool isCrystalIdxValid (unsigned i) const; + unsigned id_; + std::vector crystals_; + DiskInfo diskInfo_; + double radiusInCrystal_; + double radiusOutCrystal_; + double nominalCellSize_; + unsigned crystalOffset_; + std::shared_ptr crystalMap_; + std::vector mapToCrystal_; + std::vector crystalToMap_; + }; +} #endif diff --git a/CalorimeterGeom/inc/DiskCalorimeter.hh b/CalorimeterGeom/inc/DiskCalorimeter.hh index 5cab31c970..3b1a1e532d 100644 --- a/CalorimeterGeom/inc/DiskCalorimeter.hh +++ b/CalorimeterGeom/inc/DiskCalorimeter.hh @@ -5,68 +5,73 @@ // // Original author B. Echenard // -// Note 1: conversion of crystal <-> readout id -// readout_id = crystal_id*nRoPerCrystal ... crystal_id*nRoPerCrystal + nRoPerCrystal-1 - +#include "Offline/Mu2eInterfaces/inc/ProditionsEntity.hh" #include "Offline/CalorimeterGeom/inc/Calorimeter.hh" -#include "Offline/CalorimeterGeom/inc/CaloInfo.hh" -#include "Offline/CalorimeterGeom/inc/CaloGeomUtil.hh" +#include "Offline/CalorimeterGeom/inc/CaloG4Info.hh" #include "Offline/CalorimeterGeom/inc/Disk.hh" #include "Offline/CalorimeterGeom/inc/Crystal.hh" +#include "cetlib_except/exception.h" #include "CLHEP/Vector/ThreeVector.h" - #include #include - +#include namespace mu2e { - class DiskCalorimeter: public Calorimeter { - using DiskPtr = std::shared_ptr; - using DiskPtrs = std::vector; - using Crystals = std::vector; - using CrystalPtrs = std::vector; + class DiskCalorimeter: public Calorimeter, public ProditionsEntity { - friend class DiskCalorimeterMaker; + friend class DiskCalorimeterMaker; + public: + using Hep3Vector = CLHEP::Hep3Vector; - public: - DiskCalorimeter(); + DiskCalorimeter(); + ~DiskCalorimeter() = default; + DiskCalorimeter(const DiskCalorimeter& rhs); + DiskCalorimeter(DiskCalorimeter&& rhs) noexcept ; + DiskCalorimeter& operator=(const DiskCalorimeter& rhs) = delete; + DiskCalorimeter& operator=(DiskCalorimeter&& rhs) = delete; - size_t nDisks() const override {return disks_.size();} - const Disk& disk(size_t i) const override {return *disks_.at(i);} - const DiskPtrs& diskPtrs() const override {return disks_;} + size_t nDisks() const override {return disks_.size();} + const Disk& disk(unsigned i) const override {return disks_.at(i);} + const Disks& disks() const override {return disks_;} - size_t nCrystals() const override {return fullCrystalList_.size();} - const Crystal& crystal(size_t i) const override {return *fullCrystalList_.at(i);} - const CrystalPtrs& crystalPtrs() const override {return fullCrystalList_;} + size_t nCrystals() const override {return crystals_.size();} + const Crystal& crystal(unsigned i) const override {return *crystals_.at(i);} + const Crystals& crystals() const override {return crystals_;} - const CaloInfo& caloInfo() const override {return caloInfo_;} - const CaloGeomUtil& geomUtil() const override {return geomUtil_;} + const CaloG4Info& G4Info() const override {return G4Info_;} + void print(std::ostream &os = std::cout) const override; - const std::vector& neighbors(int crystalId) const override; - const std::vector& nextNeighbors(int crystalId) const override; - std::vector neighborsByLevel(int crystalId, int level) const override; - int crystalIdxFromPosition(const CLHEP::Hep3Vector& pos) const override; - int nearestIdxFromPosition(const CLHEP::Hep3Vector& pos) const override; + std::vector neighbors(unsigned globalId, unsigned level) const override; + bool isInsideAnyCrystal(const Hep3Vector& pos) const override; + bool isInsideAnyDisk (const Hep3Vector& pos) const override; + bool isInsideSameDisk (const Hep3Vector& front, + const Hep3Vector& back) const override; - void print(std::ostream &os = std::cout) const override; + Hep3Vector mu2eToDisk (unsigned diskId, const Hep3Vector& p) const override; + Hep3Vector diskToMu2e (unsigned diskId, const Hep3Vector& p) const override; + Hep3Vector mu2eToDiskFF (unsigned diskId, const Hep3Vector& p) const override; + Hep3Vector diskFFToMu2e (unsigned diskId, const Hep3Vector& p) const override; + Hep3Vector mu2eToCrystal(unsigned crystalId, const Hep3Vector& p) const override; + Hep3Vector crystalToMu2e(unsigned crystalId, const Hep3Vector& p) const override; + // These should REALLY go to the tracker + Hep3Vector mu2eToTracker(const Hep3Vector& p) const override; + Hep3Vector trackerToMu2e(const Hep3Vector& p) const override; - private: - double deltaZ (const CLHEP::Hep3Vector& p1, const CLHEP::Hep3Vector& p2) const; - double deltaPerp (int ic, const CLHEP::Hep3Vector& pos) const; + private: + void rebuildCrystalPtrs(); - std::vector disks_; - std::vector fullCrystalList_; //non-owning crystal pointers - CaloInfo caloInfo_; - CaloGeomUtil geomUtil_; - }; + Disks disks_; + Crystals crystals_; //non-owning cache for convenience + CaloG4Info G4Info_; + CLHEP::Hep3Vector trackerCenter_; + }; } - #endif diff --git a/CalorimeterGeom/inc/DiskGeomInfo.hh b/CalorimeterGeom/inc/DiskGeomInfo.hh deleted file mode 100644 index e694ed03d8..0000000000 --- a/CalorimeterGeom/inc/DiskGeomInfo.hh +++ /dev/null @@ -1,80 +0,0 @@ -// -// Contains gometry info of disks -// -// Original author B. Echenard -// - -#ifndef CalorimeterGeom_DiskGeomInfo_hh -#define CalorimeterGeom_DiskGeomInfo_hh - -#include "CLHEP/Vector/Rotation.h" -#include "CLHEP/Vector/ThreeVector.h" -#include - -namespace mu2e { - - class DiskGeomInfo { - - public: - DiskGeomInfo() : - size_ (CLHEP::Hep3Vector(0,0,0)), - origin_ (CLHEP::Hep3Vector(0,0,0)), - originLocal_ (CLHEP::Hep3Vector(0,0,0)), - rotation_ (CLHEP::HepRotation::IDENTITY), - inverseRotation_ (CLHEP::HepRotation::IDENTITY), - originToCrystalOrigin_(CLHEP::Hep3Vector(0,0,0)), - crystalDirection_ (CLHEP::Hep3Vector(0,0,0)), - frontFaceCenter_ (CLHEP::Hep3Vector(0,0,0)), - backFaceCenter_ (CLHEP::Hep3Vector(0,0,0)), - innerEnvelope_ (0), - outerEnvelope_ (0), - FEBZOffset_ (0), - FEBZLength_ (0) - {} - - const CLHEP::Hep3Vector& size() const {return size_; } - const CLHEP::Hep3Vector& origin() const {return origin_;} - const CLHEP::Hep3Vector& originLocal() const {return originLocal_; } - const CLHEP::Hep3Vector& originToCrystalOrigin() const {return originToCrystalOrigin_;} - const CLHEP::Hep3Vector& crystalDirection() const {return crystalDirection_;} - const CLHEP::Hep3Vector& frontFaceCenter() const {return frontFaceCenter_; } - const CLHEP::Hep3Vector& backFaceCenter() const {return backFaceCenter_; } - const CLHEP::HepRotation& rotation() const {return rotation_;} - const CLHEP::HepRotation& inverseRotation() const {return inverseRotation_;} - double innerEnvelopeR() const {return innerEnvelope_;} - double outerEnvelopeR() const {return outerEnvelope_;} - double FEBZOffset() const {return FEBZOffset_;} - double FEBZLength() const {return FEBZLength_;} - - - void size (const CLHEP::Hep3Vector& size) {size_ = size;} - void origin (const CLHEP::Hep3Vector& orig) {origin_ = orig;} - void originLocal (const CLHEP::Hep3Vector& orig) {originLocal_ = orig;} - void originToCrystalOrigin(const CLHEP::Hep3Vector& vec) {originToCrystalOrigin_ = vec;} - void crystalDirection (const CLHEP::Hep3Vector& vec) {crystalDirection_ = vec;} - void frontFaceCenter (const CLHEP::Hep3Vector& pos) {frontFaceCenter_ = pos;} - void backFaceCenter (const CLHEP::Hep3Vector& pos) {backFaceCenter_ = pos;} - void rotation (const CLHEP::HepRotation& rot) {rotation_ = rot; inverseRotation_ = rot.inverse();} - void envelopeRad (double rin, double rout) {innerEnvelope_ = rin; outerEnvelope_ = rout;} - void FEBZOffset (double val) {FEBZOffset_ = val;} - void FEBZLength (double val) {FEBZLength_ = val;} - - - private: - CLHEP::Hep3Vector size_; - CLHEP::Hep3Vector origin_; - CLHEP::Hep3Vector originLocal_; - CLHEP::HepRotation rotation_; - CLHEP::HepRotation inverseRotation_; - CLHEP::Hep3Vector originToCrystalOrigin_; - CLHEP::Hep3Vector crystalDirection_; - CLHEP::Hep3Vector frontFaceCenter_; - CLHEP::Hep3Vector backFaceCenter_; - double innerEnvelope_; - double outerEnvelope_; - double FEBZOffset_; - double FEBZLength_; - }; -} - -#endif diff --git a/CalorimeterGeom/inc/DiskInfo.hh b/CalorimeterGeom/inc/DiskInfo.hh new file mode 100644 index 0000000000..3d33bc6b40 --- /dev/null +++ b/CalorimeterGeom/inc/DiskInfo.hh @@ -0,0 +1,106 @@ +#ifndef CalorimeterGeom_DiskInfo_hh +#define CalorimeterGeom_DiskInfo_hh +// +// Contains geometry info of a disk. +// +// Design: the disk POSE (origin_ + rotation_) and the body-fixed LOCAL offsets +// are the only source of truth. Every global-frame point (face centers, crystal +// direction) and inverseRotation_ are DERIVED and cached, rebuilt in recompute_() +// whenever the pose or a local offset changes. This makes rigid moves correct by +// construction: moveDisk only needs to update the pose. +// +// NOTE for the maker: set the pose (origin + rotation) BEFORE the frontFaceCenter/ +// backFaceCenter/crystalDirection setters, since those back-solve the stored local +// offset from the global value using the current pose. +// +// Original author B. Echenard +// + +#include "CLHEP/Vector/Rotation.h" +#include "CLHEP/Vector/ThreeVector.h" +#include + +namespace mu2e { + + class DiskInfo { + + public: + DiskInfo(); + + // ---- frame conversions (single source of the transform) ---- + CLHEP::Hep3Vector toGlobal(const CLHEP::Hep3Vector& local) const; + CLHEP::Hep3Vector toLocal (const CLHEP::Hep3Vector& global) const; + CLHEP::Hep3Vector toLocalFF (const CLHEP::Hep3Vector& g) const; + CLHEP::Hep3Vector toGlobalFF(const CLHEP::Hep3Vector& l) const; + + // ---- getters (derived points return the cached value) ---- + const CLHEP::Hep3Vector& size() const {return size_;} + const CLHEP::Hep3Vector& origin() const {return origin_;} + const CLHEP::Hep3Vector& originLocal() const {return originLocal_;} + const CLHEP::Hep3Vector& originToCrystalOrigin() const {return originToCrystalOrigin_;} + const CLHEP::Hep3Vector& crystalDirection() const {return crystalDirection_;} + const CLHEP::Hep3Vector& frontFaceCenter() const {return frontFaceCenter_;} + const CLHEP::Hep3Vector& backFaceCenter() const {return backFaceCenter_;} + const CLHEP::HepRotation& rotation() const {return rotation_;} + const CLHEP::HepRotation& inverseRotation() const {return inverseRotation_;} + double crystalZLength() const {return crystalZlength_;} + double innerEnvelopeR() const {return innerEnvelope_;} + double outerEnvelopeR() const {return outerEnvelope_;} + double FEBZOffset() const {return FEBZOffset_;} + double FEBZLength() const {return FEBZLength_;} + + // ---- pose: the only mutable geometric state ---- + void origin (const CLHEP::Hep3Vector& o) {origin_ = o; recompute_();} + void rotation(const CLHEP::HepRotation& r) {rotation_ = r; recompute_();} + void setPose (const CLHEP::Hep3Vector& o, const CLHEP::HepRotation& r) + {origin_ = o; rotation_ = r; recompute_();} + + // ---- body-fixed / invariant state ---- + void originLocal (const CLHEP::Hep3Vector& o) {originLocal_ = o;} + void originToCrystalOrigin(const CLHEP::Hep3Vector& v) {originToCrystalOrigin_ = v;} // local, invariant under moves + + // ---- global-valued setters: back-solve the local offset, then rebuild ---- + void frontFaceCenter (const CLHEP::Hep3Vector& g) {ffLocal_ = toLocal(g); recompute_();} + void backFaceCenter (const CLHEP::Hep3Vector& g) {bfLocal_ = toLocal(g); recompute_();} + void crystalDirection(const CLHEP::Hep3Vector& d) {dirLocal_ = rotation_*d; recompute_();} // direction: rotate only + + // ---- frame-invariant scalars / packed sizes ---- + void size (const CLHEP::Hep3Vector& s) {size_ = s;} + void crystalZlength(double v) {crystalZlength_ = v;} + void envelopeRad (double rin, double rout) {innerEnvelope_ = rin; outerEnvelope_ = rout;} + void FEBZOffset (double v) {FEBZOffset_ = v;} + void FEBZLength (double v) {FEBZLength_ = v;} + + void print(std::ostream& os = std::cout) const; + + + private: + void recompute_(); + + // pose (source of truth) + CLHEP::Hep3Vector origin_{0,0,0}; + CLHEP::Hep3Vector originLocal_{0,0,0}; + CLHEP::HepRotation rotation_{CLHEP::HepRotation::IDENTITY}; + + // body-fixed offsets (source of truth, invariant under rigid moves) + CLHEP::Hep3Vector size_{0,0,0}; + CLHEP::Hep3Vector originToCrystalOrigin_{0,0,0}; + CLHEP::Hep3Vector ffLocal_{0,0,0}; + CLHEP::Hep3Vector bfLocal_{0,0,0}; + CLHEP::Hep3Vector dirLocal_{0,0,0}; + + // derived caches (rebuilt only in recompute_) + CLHEP::HepRotation inverseRotation_{CLHEP::HepRotation::IDENTITY}; + CLHEP::Hep3Vector frontFaceCenter_{0,0,0}; + CLHEP::Hep3Vector backFaceCenter_{0,0,0}; + CLHEP::Hep3Vector crystalDirection_{0,0,0}; + + double crystalZlength_{0}; + double innerEnvelope_{0}; + double outerEnvelope_{0}; + double FEBZOffset_{0}; + double FEBZLength_{0}; + }; +} + +#endif diff --git a/CalorimeterGeom/inc/SquareMapper.hh b/CalorimeterGeom/inc/SquareMapper.hh index b76e5aa8a6..796391b49c 100644 --- a/CalorimeterGeom/inc/SquareMapper.hh +++ b/CalorimeterGeom/inc/SquareMapper.hh @@ -1,60 +1,86 @@ -#ifndef CalorimeterGeom_SqaureMapper_hh -#define CalorimeterGeom_SqaureMapper_hh - +#ifndef CalorimeterGeom_SquareMapper_hh +#define CalorimeterGeom_SquareMapper_hh +// +// Square position map generator: +// tesselate a plane with squares starting from the center of the plane +// +// original author : Bertrand Echenard (Caltech) +// +// Use basis vector, l and k, defined as +// l = right +// k = up +// +/* + + ----------------------- + | | | | + | -1 1 | 0 1 | 1 1 | + | | | | + ----------------------- + | | | | l,k coordinates + | -1 0 | 0 0 | 1 0 | + | | | | + ----------------------- + | | | | + | -1 -1 | 0 -1 | 1 -1 | + | | | | + ----------------------- + + steps : (1,0), (0,-1), (-1,0), (0,1) (clockwise from top left corner) + segment: top=0, right=1, bottom=2,left=3 + +*/ +// Tesselation algorithm: tessalate in "rings" from the center +// for each ring, start at -l,+l (top left corner), +// then go n time each step to create the ring +// +// Neighbors add (+-1,0) or (0,+-1) +// next ring of neighbours, add +(-2,2) and go around the ring,... +// #include "Offline/CalorimeterGeom/inc/CrystalMapper.hh" #include "CLHEP/Vector/TwoVector.h" #include "CLHEP/Vector/ThreeVector.h" -#include - +#include namespace mu2e { - class SquLK { - - public: - SquLK() : l_(0),k_(0) {} - SquLK(int l, int k) : l_(l),k_(k) {} + class SquLK { + public: + constexpr SquLK() = default; + constexpr SquLK(int l, int k) : l_(l),k_(k) {} - void add(const SquLK &x) {l_+=x.l_;k_+=x.k_;} + SquLK& operator+=(const SquLK& x) {l_ += x.l_; k_ += x.k_; return *this;} - int l_; - int k_; - }; + int l_; + int k_; + }; - class SquareMapper : public CrystalMapper { + class SquareMapper : public CrystalMapper { + public: + SquareMapper() = default; - public: - SquareMapper(); + unsigned nCrystalMax(unsigned maxRing) const override; + CLHEP::Hep2Vector xyFromIndex(unsigned thisIndex) const override; + unsigned indexFromXY(double x, double y) const override; + unsigned indexFromRowCol(int row, int col) const override; + int rowFromIndex(unsigned thisIndex) const override; + int colFromIndex(unsigned thisIndex) const override; + unsigned numNeighbors(unsigned level) const override; + std::vector neighbors(unsigned thisIndex,unsigned level) const override; - int nCrystalMax (int maxRing) const override; - CLHEP::Hep2Vector xyFromIndex (int thisIndex) const override; - int indexFromXY (double x, double y) const override; - int indexFromRowCol(int nRow, int nCol) const override; - int rowFromIndex (int thisIndex) const override; - int colFromIndex (int thisIndex) const override; - bool isInsideCrystal(double x, double y, - const CLHEP::Hep3Vector& pos, - const CLHEP::Hep3Vector& size) const override; - int numNeighbors (int level) const override; - std::vector neighbors (int thisIndex, int level=1) const override; - const std::vector& apexX() const override {return apexX_;} - const std::vector& apexY() const override {return apexY_;} + private: + SquLK lk(unsigned index) const; + unsigned index(const SquLK& lk) const; + unsigned ring(const SquLK& lk) const; - - private: - SquLK lk(int index) const; - int index(const SquLK& lk) const; - int ring(const SquLK& lk) const; - - std::vector step_; - std::vector apexX_; - std::vector apexY_; - }; + inline static constexpr std::array step_{{{1,0},{0,-1},{-1,0},{0,1}}}; + inline static constexpr std::array apexX_{{-0.5,0.5,0.5,-0.5,-0.5}}; + inline static constexpr std::array apexY_{{-0.5,-0.5,0.5,0.5,-0.5}}; + }; } #endif - diff --git a/CalorimeterGeom/inc/SquareShiftMapper.hh b/CalorimeterGeom/inc/SquareShiftMapper.hh index 448d509e75..0b0ded2f90 100644 --- a/CalorimeterGeom/inc/SquareShiftMapper.hh +++ b/CalorimeterGeom/inc/SquareShiftMapper.hh @@ -1,58 +1,85 @@ -#ifndef CalorimeterGeom_SqaureShiftMapper_hh -#define CalorimeterGeom_SqaureShiftMapper_hh - +#ifndef CalorimeterGeom_SquareShiftMapper_hh +#define CalorimeterGeom_SquareShiftMapper_hh +// +// Square position map generator: +// tesselate a plane with squares, every row shifted horizontaly by 0.5 square size, +// starting from the center of the plane +// +// original author : Bertrand Echenard (Caltech) +// +// Use basis vector, l and k, defined as +// l = up right +// k = down right +// +// -------------------- +// | | | +// | 0 -1 | 1 0 | +// | | | +// | | | +// ------------------------------ +// | | | | +// | -1 -1 | 0 0 | 1 1 | l,k coordinates +// | | | | +// | | | | +// ------------------------------ +// | | | +// | -1 0 | 0 1 | +// | | | +// | | | +// -------------------- +// +// steps : (1,1) (0,1) (-1,0) (-1,-1) (0,-1) (1,0) (clockwise from top left corner) +// +// Tesselation algorithm: tessalate in "rings" from the center +// for each ring, start at 0,-l (top left corner), +// then go n time each step to create the ring +// +// Neighbors add (0,-1) and go around the ring +// next ring of neighbours, add (0,-2) and go around the ring,... +// #include "Offline/CalorimeterGeom/inc/CrystalMapper.hh" #include "CLHEP/Vector/TwoVector.h" #include "CLHEP/Vector/ThreeVector.h" -#include - +#include namespace mu2e { - class SquShiftLK { - - public: - SquShiftLK() : l_(0),k_(0) {} - SquShiftLK(int l, int k) : l_(l),k_(k) {} - - void add(const SquShiftLK &x) {l_+=x.l_;k_+=x.k_;} - - int l_; - int k_; - }; + class SquShiftLK { + public: + constexpr SquShiftLK() = default; + constexpr SquShiftLK(int l, int k) : l_(l),k_(k) {} + SquShiftLK& operator+=(const SquShiftLK& x) {l_ += x.l_; k_ += x.k_; return *this;} - class SquareShiftMapper : public CrystalMapper { + int l_; + int k_; + }; - public: - SquareShiftMapper(); - int nCrystalMax (int maxRing) const override; - CLHEP::Hep2Vector xyFromIndex (int thisIndex) const override; - int indexFromXY (double x, double y) const override; - int indexFromRowCol(int nRow, int nCol) const override; - int rowFromIndex (int thisIndex) const override; - int colFromIndex (int thisIndex) const override; - bool isInsideCrystal(double x, double y, - const CLHEP::Hep3Vector& pos, - const CLHEP::Hep3Vector& size) const override; - int numNeighbors (int level) const override; - std::vector neighbors (int thisIndex, int level=1) const override; + class SquareShiftMapper : public CrystalMapper { + public: + SquareShiftMapper() = default; - const std::vector& apexX() const override {return apexX_;} - const std::vector& apexY() const override {return apexY_;} + unsigned nCrystalMax(unsigned maxRing) const override; + CLHEP::Hep2Vector xyFromIndex(unsigned thisIndex) const override; + unsigned indexFromXY(double x, double y) const override; + unsigned indexFromRowCol(int row, int col) const override; + int rowFromIndex(unsigned thisIndex) const override; + int colFromIndex(unsigned thisIndex) const override; + unsigned numNeighbors(unsigned level) const override; + std::vector neighbors(unsigned thisIndex,unsigned level) const override; - private: - SquShiftLK lk(int index) const; - int index(const SquShiftLK& lk) const; - int ring(const SquShiftLK& lk) const; + private: + SquShiftLK lk(unsigned index) const; + unsigned index(const SquShiftLK& lk) const; + unsigned ring(const SquShiftLK& lk) const; - std::vector step_; - std::vector apexX_; - std::vector apexY_; - }; + inline static constexpr std::array step_{{{1,1},{0,1},{-1,0},{-1,-1},{0,-1},{1,0}}}; + inline static constexpr std::array apexX_{{-0.5,0.5,0.5,-0.5,-0.5}}; + inline static constexpr std::array apexY_{{-0.5,-0.5,0.5,0.5,-0.5}}; + }; } #endif diff --git a/CalorimeterGeom/src/CaloGeomUtil.cc b/CalorimeterGeom/src/CaloGeomUtil.cc deleted file mode 100644 index 690f07cf36..0000000000 --- a/CalorimeterGeom/src/CaloGeomUtil.cc +++ /dev/null @@ -1,122 +0,0 @@ -// -// Geometry and identifier info about the Calorimeter. -// -// Original author B. Echenard -// - - -#include "Offline/CalorimeterGeom/inc/CaloGeomUtil.hh" -#include "Offline/CalorimeterGeom/inc/Disk.hh" -#include "Offline/CalorimeterGeom/inc/Crystal.hh" - -#include "CLHEP/Vector/Rotation.h" -#include "CLHEP/Vector/ThreeVector.h" - - -namespace mu2e { - - - //----------------------------------------------------------------------------- - CaloGeomUtil::CaloGeomUtil(const std::vector>& disks, - const std::vector& fullCrystalList) : - disks_(disks), - fullCrystalList_(fullCrystalList), - origin_(), - trackerCenter_(), - crystalZLength_(0.0) - {} - - - //----------------------------------------------------------------------------- - CLHEP::Hep3Vector CaloGeomUtil::mu2eToCrystal(int crystalId, const CLHEP::Hep3Vector& pos) const - { - const Disk& thisDisk = disk( fullCrystalList_.at(crystalId)->diskID()); - CLHEP::Hep3Vector crysLocalPos = fullCrystalList_.at(crystalId)->localPosition(); - return thisDisk.geomInfo().rotation()*(pos-thisDisk.geomInfo().origin())-crysLocalPos; - } - - //----------------------------------------------------------------------------- - CLHEP::Hep3Vector CaloGeomUtil::mu2eToDisk(int diskId, const CLHEP::Hep3Vector& pos) const - { - const Disk& thisDisk = disk(diskId); - return (thisDisk.geomInfo().rotation())*(pos-thisDisk.geomInfo().origin()); - } - - //----------------------------------------------------------------------------- - CLHEP::Hep3Vector CaloGeomUtil::mu2eToDiskFF(int diskId, const CLHEP::Hep3Vector& pos) const - { - const Disk& thisDisk = disk(diskId); - return (thisDisk.geomInfo().rotation())*(pos-thisDisk.geomInfo().origin()) - thisDisk.geomInfo().originToCrystalOrigin(); - } - - //----------------------------------------------------------------------------- - CLHEP::Hep3Vector CaloGeomUtil::mu2eToTracker(const CLHEP::Hep3Vector& pos) const - { - return pos - trackerCenter_; - } - - //----------------------------------------------------------------------------- - CLHEP::Hep3Vector CaloGeomUtil::crystalToMu2e(int crystalId, const CLHEP::Hep3Vector& pos) const - { - const Disk& thisDisk = disk( fullCrystalList_.at(crystalId)->diskID() ); - CLHEP::Hep3Vector crysLocalPos = fullCrystalList_.at(crystalId)->localPosition(); - return thisDisk.geomInfo().inverseRotation()*(pos+crysLocalPos) + thisDisk.geomInfo().origin(); - } - - //----------------------------------------------------------------------------- - CLHEP::Hep3Vector CaloGeomUtil::diskToMu2e(int diskId, const CLHEP::Hep3Vector& pos) const - { - const Disk& thisDisk = disk(diskId); - return thisDisk.geomInfo().inverseRotation()*pos + thisDisk.geomInfo().origin(); - } - - //----------------------------------------------------------------------------- - CLHEP::Hep3Vector CaloGeomUtil::diskFFToMu2e(int diskId, const CLHEP::Hep3Vector& pos) const - { - const Disk& thisDisk = disk(diskId); - return thisDisk.geomInfo().inverseRotation()*(pos + thisDisk.geomInfo().originToCrystalOrigin()) - + thisDisk.geomInfo().origin(); - } - - //----------------------------------------------------------------------------- - CLHEP::Hep3Vector CaloGeomUtil::trackerToMu2e(const CLHEP::Hep3Vector& pos) const - { - return pos + trackerCenter_; - } - - - //----------------------------------------------------------------------------- - bool CaloGeomUtil::isInsideCalorimeter(const CLHEP::Hep3Vector& pos) const - { - for (size_t idisk=0;idisk crystalZLength_+tolerance ) return false; - if ( disk(idisk).idxFromPosition(posInSection.x(),posInSection.y()) <0) return false; - - return true; - } - - //----------------------------------------------------------------------------- - bool CaloGeomUtil::isContainedSection(const CLHEP::Hep3Vector& front, const CLHEP::Hep3Vector& back) const - { - //for (size_t idisk=0;idisk #include +#include +#include namespace mu2e { - - Disk::Disk(int id, double rCrystalIn, double rCrystalOut, double nominalCellSize, - double nominalCellLength, int offset, const CLHEP::Hep3Vector& diskOriginToCrystalOrigin) : - id_(id), - crystalList_(), - geomInfo_(), - radiusInCrystal_(rCrystalIn), - radiusOutCrystal_(rCrystalOut), - nominalCellSize_(nominalCellSize), - globalCrystalOffset_(offset), - mapToCrystal_(), - crystalToMap_(), - rowMax_(0) - { - crystalMap_ = std::shared_ptr(new SquareShiftMapper()); - - fillCrystalsIdeal(diskOriginToCrystalOrigin, nominalCellLength); //See note in DiskCalorimeterMaker - //fillCrystals(diskOriginToCrystalOrigin); //See note in DiskCalorimeterMaker - } - - - - //----------------------------------------------------------------------------- - // take the crystals from the ideal map, keeping only the crystals inside the annulus - void Disk::fillCrystalsIdeal(const CLHEP::Hep3Vector& crystalOriginInDisk, double nominalCellLength) - { - int nCrystal(0); - int nRingsMax = int(1.5*radiusOutCrystal_/nominalCellSize_); - int nCrystalMap = crystalMap_->nCrystalMax(nRingsMax); - - mapToCrystal_.insert(mapToCrystal_.begin(), nCrystalMap, invalidID_); - - for (int mapIdx=0;mapIdxxyFromIndex(mapIdx); - if (!isInsideDisk(xy.x(),xy.y(),nominalCellSize_,nominalCellSize_)) continue; - - //these crystals have been manually removed from the map by hand.... - if (std::abs(xy.x()-257.25) < 1.0 && std::abs(xy.y()-583.1) < 1.0) continue; - if (std::abs(xy.x()-257.25) < 1.0 && std::abs(xy.y()+583.1) < 1.0) continue; - if (std::abs(xy.x()+257.25) < 1.0 && std::abs(xy.y()-583.1) < 1.0) continue; - if (std::abs(xy.x()+257.25) < 1.0 && std::abs(xy.y()+583.1) < 1.0) continue; - - CLHEP::Hep3Vector size(nominalCellSize_,nominalCellSize_,nominalCellLength); - CLHEP::Hep3Vector posFF(xy.x(),xy.y(),0); - CLHEP::Hep3Vector pos = posFF + crystalOriginInDisk; - - rowMax_ = std::max(rowMax_,std::abs(crystalMap_->rowFromIndex(mapIdx))); - - mapToCrystal_[mapIdx] = nCrystal; - crystalToMap_.push_back(mapIdx); - crystalList_.push_back(Crystal(nCrystal, id_, pos, pos, size)); - - ++nCrystal; - } - } - - //----------------------------------------------------------------------------- - // fill the crystal from the real position with the measured dimensions - void Disk::fillCrystals(const CLHEP::Hep3Vector& crystalOriginInDisk) - { - /* - int nRingsMax = int(1.5*radiusOutCrystal_/nominalCellSize_); - int nCrystalMap = crystalMap_->nCrystalMax(nRingsMax); - mapToCrystal_.insert(mapToCrystal_.begin(), nCrystalMap, invalidID_); - - //yes, this will be correctly done later - WhateverDatabaseReader reader; // <-- this is the link to the database - - for (int i=0;iindexFromXY(pos.x()/nominalCellSize_,pos.y()/nominalCellSize_); - CLHEP::Hep2Vector idealPosition = nominalCellSize_*crystalMap_->xyFromIndex(mapIdx); - - mapToCrystal_[mapIdx] = i; - crystalToMap_.push_back(mapIdx); - crystalList_.push_back( Crystal(nCrystal, id_, pos, idealPosition, size) ); - } - - fixCrystalPosition(); - */ - } - - - //----------------------------------------------------------------------------- - bool Disk::isInsideDisk(double x, double y, double widthX, double widthY) const - { - std::vector apexX = crystalMap_->apexX(); - std::vector apexY = crystalMap_->apexY(); - - for (size_t i=1;i 1.0) mindist = p1.mag(); - - //farthest distance is always at an apex - double maxdist = std::max(p0.mag(),p1.mag()); - - if (mindist < radiusInCrystal_ || maxdist > radiusOutCrystal_) return false; - } - - return true; - } - - //----------------------------------------------------------------------------- - //check crystal overlaps (row after row) and move crystal around to fix boundary crossings - void Disk::fixCrystalPosition() - { - std::map> rowToCrystalId; - for (size_t cryId=0;cryIdrowFromIndex(crystalToMap_[cryId]); - rowToCrystalId[irow].push_back(cryId); - } - - double ymaxPrevious(-1e6); - for (auto& kv : rowToCrystalId){ - auto& cryList = kv.second; - - //order the crystal ids by increasing x coordinate in a given row - auto sortFunctor = [this](const auto i, const auto j) - {return crystalList_[i].localPosition().x() < - crystalList_[j].localPosition().x();}; - std::sort(std::begin(cryList),std::end(cryList),sortFunctor); - - //check if the crystals overlap in Y with the row below and fix the position - double ymax(-1e6); - for (const auto& cryIdx : cryList) { - double crYmin = crystalList_[cryIdx].localPosition().y()-crystalList_[cryIdx].size().y()/2.0; - double tolerance = crYmin - ymaxPrevious; - if (tolerance < 0) { - auto newPosition = crystalList_[cryIdx].localPosition() - CLHEP::Hep3Vector(0,tolerance,0); - crystalList_[cryIdx].setLocalPosition(newPosition); - } - ymax = std::max(ymax,crystalList_[cryIdx].localPosition().y()+crystalList_[cryIdx].size().y()/2.0); - } - ymaxPrevious = ymax; - - //check if the crystal overlap in X, and move them to the right if needed - //double origXend = crystalList_[cryList.back()].localPosition().x(); - for (size_t i=1;i0){ - // for (const auto& cryIdx : cryList) { - // auto newPosition = crystalList_[cryIdx].localPosition() - CLHEP::Hep3Vector(shiftPerCrystal,0,0); - // crystalList_[cryIdx].setLocalPosition(newPosition); - // } - //} - - // check that all crystals are still in the disk enveloppe - for (const auto& cryIdx : cryList) { - auto position = crystalList_[cryIdx].localPosition(); - auto size = crystalList_[cryIdx].size(); - if (!isInsideDisk(position.x(),position.y(),size.x(),size.y())) - throw cet::exception("Disk") << " The crystal at position="<(); + fillCrystals(crystalOrigInDisk, nominalCellLength, crystalFileName ); + } + + + //----------------------------------------------------------------------------- + // take the crystals from their measured location and correlate them with ideal map location + void Disk::fillCrystals(const CLHEP::Hep3Vector& crystalOrigInDisk, + double nominalCellLength, + const std::string& fileName) + { + //Start by reading the actual location of the crystals from the file + ConfigFileLookupPolicy configFile; + std::string fullFileName = configFile(fileName); + + std::ifstream crysFile; + crysFile.open(fullFileName); + if (!crysFile.is_open()) { + throw cet::exception("DISK_OPEN_FAILED")<<"failed to open file " << fullFileName << "\n"; + } + + std::vector realPos; + float x,y; + unsigned did,cidx; + while (crysFile >> did >> cidx >> x >> y){ + if (did != id_) continue; + if (realPos.size() != cidx) throw cet::exception("DISK_READ_FAILED")<<" invalid crystal index " << cidx << "\n"; + realPos.emplace_back(x,y,0); + } + crysFile.close(); + + + // Match the crystals to the ideal location map to facilitate navigation + unsigned nCrystal(0); + unsigned nRingsMax = unsigned(1.5*radiusOutCrystal_/nominalCellSize_); + unsigned nCrystalMap = crystalMap_->nCrystalMax(nRingsMax); + mapToCrystal_.assign(nCrystalMap, Crystal::invalidID_); + + for (unsigned mapIdx=0;mapIdxxyFromIndex(mapIdx); + const auto match = std::ranges::find_if(realPos, [&](const auto& p){ + return std::abs(p.x()-posIdeal.x()) < halfCell && std::abs(p.y()-posIdeal.y()) < halfCell; }); + if (match == realPos.end()) continue; + + CLHEP::Hep3Vector size(nominalCellSize_,nominalCellSize_,nominalCellLength); + CLHEP::Hep3Vector pos = *match + crystalOrigInDisk; + checkPosition(pos,size,nCrystal); + + mapToCrystal_[mapIdx] = nCrystal; + crystalToMap_.emplace_back(mapIdx); + crystals_.emplace_back(Crystal(nCrystal, id_, pos, size)); + ++nCrystal; + } + } + + + //----------------------------------------------------------------------------- + void Disk::checkPosition(const CLHEP::Hep3Vector& pos, const CLHEP::Hep3Vector& size, unsigned crystalID) const + { + const float overlapFraction(0.499); + for (const auto& other : crystals_) { + float dx = std::abs(other.localPosition().x()-pos.x()); + float dy = std::abs(other.localPosition().y()-pos.y()); + float sx = overlapFraction*(other.size().x() + size.x()); + float sy = overlapFraction*(other.size().y() + size.y()); + + if (dx < sx && dy < sy) throw cet::exception("DISK_POS_FAILED")<< + "crystals" << crystalID << " and " << other.localID()<<" overlap\n"; + } + } + + + //----------------------------------------------------------------------------- + std::vector Disk::idxFromRow(int thisRow) const + { + std::vector cryList; + for (size_t i=0;irowFromIndex(crystalToMap_[i]); + if (irow == thisRow) cryList.emplace_back(i); + } + return cryList; + } + + + //----------------------------------------------------------------------------- + // Position needs to be in the Disk frame + bool Disk::isInsideDisk(const CLHEP::Hep3Vector& posInDisk) const + { + if (posInDisk.z() < -diskInfo_.size().z()/2.0 || posInDisk.z() > diskInfo_.size().z()/2.0) return false; + float radius = posInDisk.perp(); + if (radius < diskInfo_.size().x() || radius > diskInfo_.size().y()) return false; + return true; + } + + + + //----------------------------------------------------------------------------- + // Position needs to be in the DiskFF frame + bool Disk::isInsideCrystal(const CLHEP::Hep3Vector& pos) const + { + const double tolerance(1e-6); + if (pos.z() < -tolerance) return false; + if (pos.z() > diskInfo_.crystalZLength()+tolerance) return false; + if (idxFromPosition(pos) == Crystal::invalidID_) return false; + return true; + } + + //----------------------------------------------------------------------------- + // Position needs to be in the DiskFF frame + unsigned Disk::idxFromPosition(const CLHEP::Hep3Vector& pos) const + { + // First, filter out obvious misses + const double tolerance(1e-6); + float perp = pos.perp(); + if (perp < diskInfo_.size().x() || + perp > diskInfo_.size().y() || + pos.z() < -tolerance || + pos.z() > diskInfo_.crystalZLength()+tolerance) return Crystal::invalidID_; + + // now look at crystal closest to ideal position - should work most of the time + auto mapIdx = crystalMap_->indexFromXY(pos.x()/nominalCellSize_,pos.y()/nominalCellSize_); + if (isCrystalIdxValid(mapIdx)) { + if (isInsideCrystal(mapToCrystal_[mapIdx],pos)) return mapToCrystal_[mapIdx]; + + // if no match, look at neighbors (likely there) + const int level(1); + const auto neighbors(crystalMap_->neighbors(mapIdx,level)); + for (const auto& idx : neighbors) { + if (isCrystalIdxValid(idx) && isInsideCrystal(mapToCrystal_[idx],pos)) return mapToCrystal_[idx]; } - - - - - - - - //----------------------------------------------------------------------------- - //calculate the bounding box surrounding a row of crystals - void Disk::boundingBoxes(int thisRow, std::vector& params) const - { - params.clear(); - - std::vector cryList; - for (size_t i=0;irowFromIndex(crystalToMap_[i]); - if (irow == thisRow) cryList.push_back(i); - } - if (cryList.empty()) return; - - auto sortFunctor = [this](const auto i, const auto j) - {return crystalList_[i].localPosition().x() < - crystalList_[j].localPosition().x();}; - std::sort(std::begin(cryList),std::end(cryList),sortFunctor); - - double ymin(1e6),ymax(-1e6),xInmin(-1e6),xImax(1e6); - for (size_t i=0;i 2*crystalList_[cryList[i]].size().x()){xInmin = x0; xImax=x1;} - } - - double xOutmin = crystalList_[cryList.front()].localPosition().x()-crystalList_[cryList.front()].size().x()/2.0; - double xOutmax = crystalList_[cryList.back()].localPosition().x() +crystalList_[cryList.back()].size().x()/2.0; - params.push_back(xOutmin); - params.push_back(xOutmax); - params.push_back(ymin); - params.push_back(ymax); - if (xInmin>-1e5) {params.push_back(xInmin);params.push_back(xImax);} - } - - - //----------------------------------------------------------------------------- - int Disk::idMinCrystalInside(int row) const - { - if (std::abs(row) > rowMax_) return invalidID_; - int idx(0); - while (mapToCrystal_[crystalMap_->indexFromRowCol(row,idx)] == invalidID_) ++idx; - return mapToCrystal_[crystalMap_->indexFromRowCol(row,idx)]; - } - - //----------------------------------------------------------------------------- - int Disk::idMaxCrystalInside(int row) const - { - if (std::abs(row) > rowMax_) return invalidID_; - int idx = int(radiusOutCrystal_/nominalCellSize_)+2; - while ( !isMapIdxValid(crystalMap_->indexFromRowCol(row,idx)) ) --idx; - while ( mapToCrystal_[crystalMap_->indexFromRowCol(row,idx)]==invalidID_ && idx>0) --idx; - return mapToCrystal_[crystalMap_->indexFromRowCol(row,idx)]; - } - - - - //----------------------------------------------------------------------------- - int Disk::idxFromPosition(double x, double y) const - { - // this only works for ideal crystal placement - //int mapIdx = crystalMap_->indexFromXY(x/nominalCellSize_,y/nominalCellSize_); - //if (mapIdx < mapToCrystal_.size() ) return mapToCrystal_.at(mapIdx); - //return invalidID_; - - int mapIdx = crystalMap_->indexFromXY(x/nominalCellSize_,y/nominalCellSize_); - - if (isCrystalIdxValid(mapIdx) && isInsideCrystal(mapToCrystal_[mapIdx],x,y)) return mapToCrystal_[mapIdx]; - - const int level(1); - const auto neighbors(crystalMap_->neighbors(mapIdx,level)); - for (const auto& idx : neighbors) { - if (isCrystalIdxValid(idx) && isInsideCrystal(mapToCrystal_[idx],x,y)) return mapToCrystal_[idx]; - } - - return invalidID_; - } - - - - - //----------------------------------------------------------------------------- - bool Disk::isInsideCrystal(int icry, double x, double y) const - { - if (icry >= static_cast(crystalList_.size())) return false; - return crystalMap_->isInsideCrystal(x, y, crystalList_[icry].localPosition(), crystalList_[icry].size()); - } - - - //----------------------------------------------------------------------------- - //find the local indexes of the crystal neighbors for a given level (level = number of rings away) - std::vector Disk::findLocalNeighbors(int crystalId, int level) const - { - std::vector list; - std::vector temp(crystalMap_->neighbors(crystalToMap_.at(crystalId),level)); - - for (const auto& mapIdx : temp) { - if (isCrystalIdxValid(mapIdx)) list.push_back(mapToCrystal_.at(mapIdx)); - } - return list; - } - - - //----------------------------------------------------------------------------- - //find the nearest crystals from the position based on ideal mapping - std::vector Disk::nearestNeighborsFromPos(double x, double y) const - { - int level(1); - std::vector list; - - auto mapIdx = crystalMap_->indexFromXY(x/nominalCellSize_,y/nominalCellSize_); - if (isCrystalIdxValid(mapIdx)) list.push_back(mapToCrystal_.at(mapIdx)); - - while (list.empty()) { - for (auto mapIdx : crystalMap_->neighbors(mapIdx,level)) { - if (isCrystalIdxValid(mapIdx)) list.push_back(mapToCrystal_[mapIdx]); - } - ++level; - } - return list; - } - - - - //----------------------------------------------------------------------------- - const bool Disk::isCrystalIdxValid(int i) const { - return i < static_cast(mapToCrystal_.size()) && mapToCrystal_[i]!=invalidID_; - } - - //----------------------------------------------------------------------------- - const bool Disk::isMapIdxValid(int i) const { - return i < static_cast(mapToCrystal_.size()); - } - - - //----------------------------------------------------------------------------- - void Disk::print(std::ostream &os) const - { - os<<"Disk "<= crystals_.size()) return false; + + const auto& cpos = crystals_[icry].localPosition(); + const auto& csize = crystals_[icry].size(); + if (std::abs(cpos.x()-pos.x()) > 0.5*csize.x() ) return false; + if (std::abs(cpos.y()-pos.y()) > 0.5*csize.y() ) return false; + if (pos.z() < -1e-6 || pos.z() > csize.z() ) return false; + + return true; + } + + //----------------------------------------------------------------------------- + //find the indices of the crystal neighbors for a given level (level = number of rings away) + // local in, local out + std::vector Disk::neighbors(unsigned localId, unsigned level) const + { + std::vector list; + const auto mapNbrs = crystalMap_->neighbors(crystalToMap_.at(localId), level); + for (const auto& mapIdx : mapNbrs) + if (isCrystalIdxValid(mapIdx)) list.emplace_back(mapToCrystal_.at(mapIdx)); + return list; + } + + //----------------------------------------------------------------------------- + bool Disk::isCrystalIdxValid(unsigned i) const { + return i < mapToCrystal_.size() && mapToCrystal_[i] != Crystal::invalidID_; + } + + //----------------------------------------------------------------------------- + //Move and rotate the disk. The displacement and rotation are relatives + void Disk::moveDisk(const CLHEP::Hep3Vector& shift, const CLHEP::HepRotation& rot) + { + diskInfo_.setPose(diskInfo_.origin() + shift, diskInfo_.rotation() * rot); + for (auto& c : crystals_) c.setPosition(diskInfo_.toGlobal(c.localPosition())); + } + + //----------------------------------------------------------------------------- + void Disk::moveCrystal(unsigned localId, const CLHEP::Hep3Vector& shift) + { + auto& crystal = crystals_.at(localId); + crystal.setLocalPosition(crystal.localPosition() + shift); + crystal.setPosition(diskInfo_.toGlobal(crystal.localPosition())); + } + + + //----------------------------------------------------------------------------- + CLHEP::Hep3Vector Disk::mu2eToDisk (const CLHEP::Hep3Vector& p) const { + return diskInfo_.toLocal(p); + } + + CLHEP::Hep3Vector Disk::diskToMu2e (const CLHEP::Hep3Vector& p) const { + return diskInfo_.toGlobal(p); + } + + CLHEP::Hep3Vector Disk::mu2eToDiskFF(const CLHEP::Hep3Vector& p) const { + return diskInfo_.toLocalFF(p); + } + + CLHEP::Hep3Vector Disk::diskFFToMu2e(const CLHEP::Hep3Vector& p) const { + return diskInfo_.toGlobalFF(p); + } + + CLHEP::Hep3Vector Disk::mu2eToCrystal(unsigned localId, const CLHEP::Hep3Vector& p) const { + return diskInfo_.toLocal(p) - crystal(localId).localPosition(); + } + + CLHEP::Hep3Vector Disk::crystalToMu2e(unsigned localId, const CLHEP::Hep3Vector& p) const { + return diskInfo_.toGlobal(p + crystal(localId).localPosition()); + } + + + //----------------------------------------------------------------------------- + void Disk::print(std::ostream &os) const + { + os<<"Disk "<& DiskCalorimeter::neighbors(int crystalId) const - { - return fullCrystalList_.at(crystalId)->neighbors(); - } - - - //----------------------------------------------------------------------------- - const std::vector& DiskCalorimeter::nextNeighbors(int crystalId) const - { - return fullCrystalList_.at(crystalId)->nextNeighbors(); - } - - - //----------------------------------------------------------------------------- - std::vector DiskCalorimeter::neighborsByLevel(int crystalId, int level) const - { - int iv = fullCrystalList_.at(crystalId)->diskID(); - int offset = disk(iv).crystalOffset(); - - std::vector list = disk(iv).findLocalNeighbors(fullCrystalList_.at(crystalId)->localID(),level); - transform(list.begin(), list.end(), list.begin(),[=](int i){return (i<0) ? i : i+offset;}); - - return list; - } - - - //----------------------------------------------------------------------------- - int DiskCalorimeter::crystalIdxFromPosition(const CLHEP::Hep3Vector& pos) const - { - for (size_t idisk=0;idiskdeltaZ(disk1->geomInfo().origin(),pos) < - this->deltaZ(disk2->geomInfo().origin(),pos);}; - auto minDistIter = std::min_element(disks_.begin(),disks_.end(),minPred); - - std::vector cand = disks_[0]->nearestNeighborsFromPos(posInSection.x(),posInSection.y()); - auto bestPred = [&](int ic1, int ic2) {return this->deltaPerp(ic1,pos) < this->deltaPerp(ic2,pos);}; - auto bestCandIter = std::min_element(cand.begin(),cand.end(),bestPred); - - return *bestCandIter+(*minDistIter)->crystalOffset(); + } + + std::vector DiskCalorimeter::neighbors(unsigned globalId, unsigned level) const + { + const Crystal& c = crystal(globalId); + const Disk& d = disk(c.diskID()); + auto ids = d.neighbors(c.localID(), level); + for (auto& v : ids) v += d.crystalOffset(); + return ids; + } + + bool DiskCalorimeter::isInsideAnyDisk(const CLHEP::Hep3Vector& pos) const + { + for (const auto& disk : disks_) { + CLHEP::Hep3Vector posInDisk = mu2eToDisk(disk.id(),pos); + if (disk.isInsideDisk(posInDisk)) return true; } - - - //----------------------------------------------------------------------------- - double DiskCalorimeter::deltaZ(const CLHEP::Hep3Vector& p1, const CLHEP::Hep3Vector& p2) const - { - return std::abs(p1.z()-p2.z()); - } - - - //----------------------------------------------------------------------------- - double DiskCalorimeter::deltaPerp(int ic, const CLHEP::Hep3Vector& pos) const - { - return sqrt((fullCrystalList_.at(ic)->position()-pos).perp2()); + return false; + } + + bool DiskCalorimeter::isInsideAnyCrystal(const CLHEP::Hep3Vector& pos) const + { + for (const auto& disk : disks_) { + CLHEP::Hep3Vector posInDisk = mu2eToDiskFF(disk.id(),pos); + if (disk.isInsideCrystal(posInDisk)) return true; } - - - //----------------------------------------------------------------------------- - void DiskCalorimeter::print(std::ostream &os) const - { - os<<"Disk calorimeter "< +#include namespace mu2e { - - SquareMapper::SquareMapper() : - step_(), - apexX_({-0.5,0.5,0.5,-0.5,-0.5}), - apexY_({-0.5,-0.5,0.5,0.5,-0.5}) - { - step_.push_back( SquLK( 1, 0) ); //right - step_.push_back( SquLK( 0,-1) ); //down - step_.push_back( SquLK(-1, 0) ); //left - step_.push_back( SquLK( 0, 1) ); //up - } - - //-------------------------------------------------------------------------------- - int SquareMapper::nCrystalMax(int maxRing) const {return (2*maxRing+1)*(2*maxRing+1);} - - - //-------------------------------------------------------------------------------- - CLHEP::Hep2Vector SquareMapper::xyFromIndex(int thisIndex) const - { - SquLK thisLK = lk(thisIndex); - return CLHEP::Hep2Vector(thisLK.l_,thisLK.k_); - } - - - int SquareMapper::indexFromXY(double x0, double y0) const - { - int l = int(std::abs(x0)+0.5); - int k = int(std::abs(y0)+0.5); - if (x0<0) l *= -1; - if (y0<0) k *= -1; - - SquLK lk(l,k); - return index(lk); - } - - - //-------------------------------------------------------------------------------- - int SquareMapper::indexFromRowCol(int nRow, int nCol) const - { - SquLK lk(nCol,nRow); - return index(lk); - } - - int SquareMapper::rowFromIndex(int thisIndex) const - { - SquLK thisLK = lk(thisIndex); - return thisLK.k_; - } - - int SquareMapper::colFromIndex(int thisIndex) const - { - SquLK thisLK = lk(thisIndex); - return thisLK.l_; - } - - - //-------------------------------------------------------------------------------- - bool SquareMapper::isInsideCrystal(double x, double y, const CLHEP::Hep3Vector& pos, - const CLHEP::Hep3Vector& size) const - { - return (std::abs(x-pos.x()) < 0.5*size.x()) && (std::abs(y-pos.y()) < 0.5*size.y()); - } - - //-------------------------------------------------------------------------------- - std::vector SquareMapper::neighbors(int thisIndex, int level) const - { - if (level<1) return std::vector{}; - - std::vector thisNeighbour; - - SquLK init = lk(thisIndex); - SquLK lk(init.l_ - level, init.k_ + level); - - for (size_t i=0;i -nRing) pos += 2*nRing + nRing - thisLK.k_; - if ( thisLK.k_ == -nRing && thisLK.l_ > -nRing) pos += 4*nRing + nRing - thisLK.l_; - if ( thisLK.l_ == -nRing && thisLK.k_ < nRing) pos += 6*nRing + nRing + thisLK.k_; - return pos; - } - - //-------------------------------------------------------------------------------- - int SquareMapper::ring(const SquLK &thisLK) const - { - return std::max(std::abs(thisLK.l_),std::abs(thisLK.k_)); - } - - //-------------------------------------------------------------------------------- - int SquareMapper::numNeighbors(int level) const - { - return 2*level*static_cast(step_.size()); - } + //-------------------------------------------------------------------------------- + unsigned SquareMapper::nCrystalMax(unsigned maxRing) const + { + return (2u*maxRing+1u)*(2u*maxRing+1u); + } + + //-------------------------------------------------------------------------------- + CLHEP::Hep2Vector SquareMapper::xyFromIndex(unsigned thisIndex) const + { + SquLK thisLK = lk(thisIndex); + return CLHEP::Hep2Vector(thisLK.l_,thisLK.k_); + } + + //-------------------------------------------------------------------------------- + unsigned SquareMapper::indexFromXY(double x0, double y0) const + { + int l = int(std::abs(x0)+0.5); + int k = int(std::abs(y0)+0.5); + + if (x0 < 0) l = -l; + if (y0 < 0) k = -k; + + return index({l,k}); + } + + //-------------------------------------------------------------------------------- + unsigned SquareMapper::indexFromRowCol(int row, int col) const + { + return index({col,row}); + } + + //-------------------------------------------------------------------------------- + int SquareMapper::rowFromIndex(unsigned thisIndex) const + { + return lk(thisIndex).k_; + } + + //-------------------------------------------------------------------------------- + int SquareMapper::colFromIndex(unsigned thisIndex) const + { + return lk(thisIndex).l_; + } + + //-------------------------------------------------------------------------------- + std::vector + SquareMapper::neighbors(unsigned thisIndex, unsigned level) const + { + if (!level) + return {}; + + std::vector result; + result.reserve(numNeighbors(level)); + + auto current = lk(thisIndex); + current.l_ -= int(level); + current.k_ += int(level); + + for (const auto& s : step_) + for (unsigned j=0; j<2u*level; ++j) { current += s; result.push_back(index(current)); } + return result; + } + + //-------------------------------------------------------------------------------- + SquLK SquareMapper::lk(unsigned thisIndex) const + { + if (!thisIndex) return {0,0}; + + const unsigned nRing = unsigned(0.5*sqrt(double(thisIndex))+0.5); + const unsigned first = (2u*nRing-1u)*(2u*nRing-1u); + const unsigned nSeg = (thisIndex-first)/(2u*nRing); + const unsigned nPos = (thisIndex-first)%(2u*nRing); + + const int r = int(nRing); + const int p = int(nPos); + if (nSeg==0) return {-r+p, r}; + if (nSeg==1) return { r, r-p}; + if (nSeg==2) return { r-p, -r}; + return {-r, -r+p}; + } + + //-------------------------------------------------------------------------------- + unsigned SquareMapper::index(const SquLK& thisLK) const + { + if (thisLK.l_==0 && thisLK.k_==0) return 0u; + + unsigned nRing = ring(thisLK); + unsigned pos = (2u*nRing-1u)*(2u*nRing-1u); + int r = int(nRing); + + // add position along segment + if (thisLK.k_ == r && thisLK.l_ < r) pos += unsigned(r + thisLK.l_); + if (thisLK.l_ == r && thisLK.k_ > -r) pos += unsigned(3*r - thisLK.k_); + if (thisLK.k_ == -r && thisLK.l_ > -r) pos += unsigned(5*r - thisLK.l_); + if (thisLK.l_ == -r && thisLK.k_ < r) pos += unsigned(7*r + thisLK.k_); + return pos; + } + + //-------------------------------------------------------------------------------- + unsigned SquareMapper::ring(const SquLK& thisLK) const + { + return unsigned(std::max(std::abs(thisLK.l_), std::abs(thisLK.k_))); + } + + //-------------------------------------------------------------------------------- + unsigned SquareMapper::numNeighbors(unsigned level) const + { + return 2u*level*unsigned(step_.size()); + } } diff --git a/CalorimeterGeom/src/SquareShiftMapper.cc b/CalorimeterGeom/src/SquareShiftMapper.cc index 03252d3020..6c20552563 100644 --- a/CalorimeterGeom/src/SquareShiftMapper.cc +++ b/CalorimeterGeom/src/SquareShiftMapper.cc @@ -1,203 +1,135 @@ -// -// Square position map generator: -// tesselate a plane with squares, every row shifted horizontaly by 0.5 square size, -// starting from the center of the plane -// -// original author : Bertrand Echenard (Caltech) -// -// Use basis vector, l and k, defined as -// l = up right -// k = down right -// -// -------------------- -// | | | -// | 0 -1 | 1 0 | -// | | | -// | | | -// ------------------------------ -// | | | | -// | -1 -1 | 0 0 | 1 1 | l,k coordinates -// | | | | -// | | | | -// ------------------------------ -// | | | -// | -1 0 | 0 1 | -// | | | -// | | | -// -------------------- -// -// steps : (1,1) (0,1) (-1,0) (-1,-1) (0,-1) (1,0) (clockwise from top left corner) -// -// Tesselation algorithm: tessalate in "rings" from the center -// for each ring, start at 0,-l (top left corner), -// then go n time each step to create the ring -// -// Neighbors add (0,-1) and go around the ring -// next ring of neighbours, add (0,-2) and go around the ring,... -// - #include "Offline/CalorimeterGeom/inc/SquareShiftMapper.hh" #include "CLHEP/Vector/TwoVector.h" #include "CLHEP/Vector/ThreeVector.h" +#include +#include namespace mu2e { - - SquareShiftMapper::SquareShiftMapper() : - step_(), - apexX_({-0.5, 0.5,0.5,-0.5,-0.5}), - apexY_({-0.5,-0.5,0.5, 0.5,-0.5}) - { - step_.push_back( SquShiftLK( 1, 1) ); //right - step_.push_back( SquShiftLK( 0, 1) ); //down right - step_.push_back( SquShiftLK(-1, 0) ); //down left - step_.push_back( SquShiftLK(-1,-1) ); //left - step_.push_back( SquShiftLK( 0,-1) ); //up left - step_.push_back( SquShiftLK( 1, 0) ); //up right - } - - - //-------------------------------------------------------------------------------- - int SquareShiftMapper::nCrystalMax(int maxRing) const {return 3*maxRing*(maxRing+1)+1;} - - - - //-------------------------------------------------------------------------------- - CLHEP::Hep2Vector SquareShiftMapper::xyFromIndex(int thisIndex) const - { - SquShiftLK thisLK = lk(thisIndex); - return CLHEP::Hep2Vector( (thisLK.l_+thisLK.k_)/2.0, (thisLK.l_-thisLK.k_) ); - } - - int SquareShiftMapper::indexFromXY(double x0, double y0) const - { - int l(0),k(0); - int ny = (y0>0) ? int(std::abs(y0)+0.5) : -int(std::abs(y0)+0.5); - - if (ny%2==0) - { - int nx = (x0>0) ? int(std::abs(x0)+0.5) : -int(std::abs(x0)+0.5); - l = nx+ny/2; - k = nx-ny/2; - } else { - int nx = (x0>0) ? int(std::abs(x0)) : (-int(std::abs(x0))-1); - l = nx + (ny+1)/2; - k = nx - (ny-1)/2; - } - - SquShiftLK lk(l,k); - return index(lk); - } - - - - //-------------------------------------------------------------------------------- - int SquareShiftMapper::indexFromRowCol(int nRow, int nCol) const - { - int k = nRow/2-nRow + nCol; - int l = nRow/2 + nCol; - - SquShiftLK lk(l,k); - return index(lk); - } - - int SquareShiftMapper::rowFromIndex(int thisIndex) const - { - SquShiftLK thisLK = lk(thisIndex); - return thisLK.l_ - thisLK.k_; - } - - int SquareShiftMapper::colFromIndex(int thisIndex) const - { - SquShiftLK thisLK = lk(thisIndex); - if ((thisLK.l_+thisLK.k_)%2) return (thisLK.l_+thisLK.k_)/2; - else if (thisLK.l_!=0) return (thisLK.l_+thisLK.k_-abs(thisLK.l_)/thisLK.l_)/2+abs(thisLK.l_)/thisLK.l_; - else return thisLK.k_/2+abs(thisLK.k_)/thisLK.k_; - } - - - //-------------------------------------------------------------------------------- - bool SquareShiftMapper::isInsideCrystal(double x, double y, const CLHEP::Hep3Vector& pos, - const CLHEP::Hep3Vector& size) const - { - return (std::abs(x-pos.x()) < 0.5*size.x()) && (std::abs(y-pos.y()) < 0.5*size.y()); - } - - - - //-------------------------------------------------------------------------------- - std::vector SquareShiftMapper::neighbors(int thisIndex, int level) const - { - if (level<1) return std::vector{}; - - std::vector thisNeighbour; - - SquShiftLK init = lk(thisIndex); - SquShiftLK lk(init.l_, init.k_ - level); - - for (size_t i=0; i(thisIndex-1)/3.0)); - int nSeg = (thisIndex - 3*nRing*(nRing-1)-1) / nRing; - int nPos = (thisIndex - 3*nRing*(nRing-1)-1) % nRing; - - if (nSeg<0) return SquShiftLK(0,0); - int l = nPos*step_[nSeg].l_; - int k = -nRing + nPos*step_[nSeg].k_; - - for (int i=0;i thisLK.k_) return pos + thisLK.l_; - return pos + 3*nRing +std::abs(thisLK.l_); - } - - - //-------------------------------------------------------------------------------- - int SquareShiftMapper::ring(const SquShiftLK& thisLK) const - { - if (thisLK.l_*thisLK.k_>0) return std::max(std::abs(thisLK.l_),std::abs(thisLK.k_)); - return std::abs(thisLK.l_-thisLK.k_); - } - - //-------------------------------------------------------------------------------- - int SquareShiftMapper::numNeighbors(int level) const - { - return level*static_cast(step_.size()); - } + //-------------------------------------------------------------------------------- + unsigned SquareShiftMapper::nCrystalMax(unsigned maxRing) const + { + return 3u*maxRing*(maxRing+1u)+1u; + } + + //-------------------------------------------------------------------------------- + CLHEP::Hep2Vector SquareShiftMapper::xyFromIndex(unsigned thisIndex) const + { + SquShiftLK thisLK = lk(thisIndex); + return CLHEP::Hep2Vector( (thisLK.l_+thisLK.k_)/2.0, (thisLK.l_-thisLK.k_) ); + } + + //-------------------------------------------------------------------------------- + unsigned SquareShiftMapper::indexFromXY(double x0, double y0) const + { + int l=0, k=0; + int ny = (y0>0) ? int(std::abs(y0)+0.5) : -int(std::abs(y0)+0.5); + + if (ny%2==0) { + int nx = (x0>0) ? int(std::abs(x0)+0.5) : -int(std::abs(x0)+0.5); + l = nx + ny/2; + k = nx - ny/2; + } else { + int nx = (x0>0) ? int(std::abs(x0)) : -int(std::abs(x0))-1; + l = nx + (ny+1)/2; + k = nx - (ny-1)/2; + } + + return index({l,k}); + } + + unsigned SquareShiftMapper::indexFromRowCol(int row, int col) const + { + const int p = ((row % 2) + 2) % 2; // 0 if row even, 1 if odd (negatives safe) + const int s = 2*col + p; // s = l + k + const int l = (s + row)/2; + const int k = (s - row)/2; + return index({l,k}); + } + + //-------------------------------------------------------------------------------- + int SquareShiftMapper::rowFromIndex(unsigned thisIndex) const + { + auto lk0 = lk(thisIndex); + return lk0.l_-lk0.k_; + } + + //-------------------------------------------------------------------------------- + int SquareShiftMapper::colFromIndex(unsigned thisIndex) const + { + auto lk0 = lk(thisIndex); + int s = lk0.l_ + lk0.k_; + return (s >= 0) ? s/2 : -(((-s) + 1)/2); // floor((l+k)/2) + } + + //-------------------------------------------------------------------------------- + std::vector SquareShiftMapper::neighbors(unsigned thisIndex, unsigned level) const + { + if (!level) return {}; + + std::vector result; + result.reserve(numNeighbors(level)); + + auto current = lk(thisIndex); + current.k_ -= int(level); + for (const auto& s : step_) + for (unsigned j=0; j=step_.size()) return {0,0}; + + int l = int(nPos)*step_[nSeg].l_; + int k = -int(nRing)+int(nPos)*step_[nSeg].k_; + + for (unsigned i=0; ilk0.k_) return pos + unsigned(lk0.l_); + return pos + 3u*nRing + unsigned(std::abs(lk0.l_)); + } + + //-------------------------------------------------------------------------------- + unsigned SquareShiftMapper::ring(const SquShiftLK& lk0) const + { + if (lk0.l_*lk0.k_>0) + return unsigned(std::max(std::abs(lk0.l_),std::abs(lk0.k_))); + return unsigned(std::abs(lk0.l_-lk0.k_)); + } + + //-------------------------------------------------------------------------------- + unsigned SquareShiftMapper::numNeighbors(unsigned level) const + { + return level*unsigned(step_.size()); + } } diff --git a/EventDisplay/src/DataInterface.cc b/EventDisplay/src/DataInterface.cc index d2cd468664..39fb212ed0 100644 --- a/EventDisplay/src/DataInterface.cc +++ b/EventDisplay/src/DataInterface.cc @@ -360,30 +360,30 @@ void DataInterface::fillGeometry() { mu2e::GeomHandle calo; - double diskCaseDZLength = calo->caloInfo().getDouble("diskCaseZLength")/2.0; - double diskRadiusIn = calo->caloInfo().getDouble("caloDiskRadiusIn"); - double diskRadiusOut = calo->caloInfo().getDouble("caloDiskRadiusOut"); - - double FPCarbonDZ = calo->caloInfo().getDouble("FPCarbonZLength")/2.0; - double FPFoamDZ = calo->caloInfo().getDouble("FPFoamZLength")/2.0; - double FPCoolPipeRadius = calo->caloInfo().getDouble("FPCoolPipeRadius"); - double pipeRadius = calo->caloInfo().getDouble("pipeRadius"); + double diskCaseDZLength = calo->G4Info().get("diskCaseZLength")/2.0; + double diskRadiusIn = calo->G4Info().get("caloDiskRadiusIn"); + double diskRadiusOut = calo->G4Info().get("caloDiskRadiusOut"); + + double FPCarbonDZ = calo->G4Info().get("FPCarbonZLength")/2.0; + double FPFoamDZ = calo->G4Info().get("FPFoamZLength")/2.0; + double FPCoolPipeRadius = calo->G4Info().get("FPCoolPipeRadius"); + double pipeRadius = calo->G4Info().get("pipeRadius"); double frontPanelHalfThick = (2.0*FPCarbonDZ+2.0*FPFoamDZ-pipeRadius+FPCoolPipeRadius)/2.0; - double crystalDXY = calo->caloInfo().getDouble("crystalXYLength")/2.0; - double crystalDZ = calo->caloInfo().getDouble("crystalZLength")/2.0; - double crystalFrameDZ = calo->caloInfo().getDouble("crystalCapZLength")/2.0; - double wrapperHalfThick = calo->caloInfo().getDouble("wrapperThickness")/2.0; + double crystalDXY = calo->G4Info().get("crystalXYLength")/2.0; + double crystalDZ = calo->G4Info().get("crystalZLength")/2.0; + double crystalFrameDZ = calo->G4Info().get("crystalCapZLength")/2.0; + double wrapperHalfThick = calo->G4Info().get("wrapperThickness")/2.0; double wrapperDXY = crystalDXY + 2.0*wrapperHalfThick; double wrapperDZ = crystalDZ + 2.0*crystalFrameDZ; - double FEEDZ = calo->caloInfo().getDouble("FEEZLength")/2.0; - double FEEBoxThickness = calo->caloInfo().getDouble("FEEBoxThickness"); + double FEEDZ = calo->G4Info().get("FEEZLength")/2.0; + double FEEBoxThickness = calo->G4Info().get("FEEBoxThickness"); double FEEBoxDZ = FEEDZ + 2*FEEBoxThickness; - double BPPipeRadiusHigh = calo->caloInfo().getDouble("BPPipeRadiusHigh"); - double BPPipeDZOffset = calo->caloInfo().getDouble("BPPipeZOffset")/2.0; + double BPPipeRadiusHigh = calo->G4Info().get("BPPipeRadiusHigh"); + double BPPipeDZOffset = calo->G4Info().get("BPPipeZOffset")/2.0; double BPFEEDZ = FEEBoxDZ + BPPipeDZOffset + BPPipeRadiusHigh; - double holeDZ = calo->caloInfo().getDouble("BPHoleZLength")/2.0; + double holeDZ = calo->G4Info().get("BPHoleZLength")/2.0; double zHalfBP = BPFEEDZ+holeDZ; double crystalDiskLogOffset = frontPanelHalfThick - zHalfBP; @@ -391,7 +391,7 @@ void DataInterface::fillGeometry() size_t icrystal=0; for(unsigned int idisk=0; idisknDisks(); idisk++) { - CLHEP::Hep3Vector diskPos = calo->disk(idisk).geomInfo().origin() - _detSysOrigin; + CLHEP::Hep3Vector diskPos = calo->disk(idisk).diskInfo().origin() - _detSysOrigin; diskPos += CLHEP::Hep3Vector(0.0, 0.0, crystalDiskLogOffset); findBoundaryP(_calorimeterMinmax, diskPos.x()+diskRadiusOut, diskPos.y()+diskRadiusOut, diskPos.z()+diskCaseDZLength); diff --git a/EventGenerator/src/CORSIKAEventGenerator_module.cc b/EventGenerator/src/CORSIKAEventGenerator_module.cc index 4f79e95dc3..ca2702d695 100644 --- a/EventGenerator/src/CORSIKAEventGenerator_module.cc +++ b/EventGenerator/src/CORSIKAEventGenerator_module.cc @@ -77,7 +77,7 @@ namespace mu2e { std::vector _worldIntersections; - bool _geomInfoObtained = false; + bool _diskInfoObtained = false; std::string _corsikaModuleLabel = "FromCorsikaBinary"; float _worldXmin = 0; @@ -148,7 +148,7 @@ namespace mu2e { _stoff = *stoH; } - if (!_geomInfoObtained) { + if (!_diskInfoObtained) { GeomHandle env; GeomHandle worldGeom; GeomHandle detsys; @@ -177,12 +177,12 @@ namespace mu2e { { GeomHandle calorimeter; _cosmicReferencePointInMu2e = Hep3Vector(detsys->getOrigin().x(), - 0, calorimeter->disk(0).geomInfo().origin().z()); + 0, calorimeter->disk(0).diskInfo().origin().z()); } else if (_refPointChoice == "UNDEFINED") _cosmicReferencePointInMu2e = Hep3Vector(0., 0, 0.); - _geomInfoObtained = true; + _diskInfoObtained = true; } art::Handle corsikaParticles; diff --git a/EventGenerator/src/CaloCalibGun_module.cc b/EventGenerator/src/CaloCalibGun_module.cc index 9d37b8ca7f..c61c56f772 100644 --- a/EventGenerator/src/CaloCalibGun_module.cc +++ b/EventGenerator/src/CaloCalibGun_module.cc @@ -19,7 +19,7 @@ #include "Offline/DataProducts/inc/PDGCode.hh" #include "Offline/SeedService/inc/SeedService.hh" #include "Offline/MCDataProducts/inc/ProcessCode.hh" -#include "Offline/CalorimeterGeom/inc/DiskCalorimeter.hh" +#include "Offline/CalorimeterGeom/inc/Calorimeter.hh" #include "Offline/Mu2eUtilities/inc/RandomUnitSphere.hh" #include "Offline/MCDataProducts/inc/GenId.hh" #include "Offline/MCDataProducts/inc/GenParticle.hh" @@ -129,22 +129,22 @@ namespace mu2e { } void CaloCalibGun::beginRun(art::Run&){ - const DiskCalorimeter *_cal = GeomHandle().get(); + const Calorimeter *_cal = GeomHandle().get(); - _pipeRadius = _cal->caloInfo().getDouble("pipeRadius"); - _pipeTorRadius = _cal->caloInfo().getVDouble("pipeTorRadius"); - _zPipeCenter = _cal->disk(_nDisk).geomInfo().origin()-CLHEP::Hep3Vector(0,0,_cal->disk(_nDisk).geomInfo().size().z()/2.0-_pipeRadius); - _nPipes = _cal->caloInfo().getInt("nPipes"); + _pipeRadius = _cal->G4Info().get("pipeRadius"); + _pipeTorRadius = _cal->G4Info().get>("pipeTorRadius"); + _zPipeCenter = _cal->disk(_nDisk).diskInfo().origin()-CLHEP::Hep3Vector(0,0,_cal->disk(_nDisk).diskInfo().size().z()/2.0-_pipeRadius); + _nPipes = _cal->G4Info().get("nPipes"); //Define the parameters of the pipes: - phi_lbd = _cal->caloInfo().getVDouble("largeTorPhi"); - phi_sbd = _cal->caloInfo().getVDouble("smallTorPhi"); - phi_end = _cal->caloInfo().getVDouble("straightEndPhi"); - ysmall = _cal->caloInfo().getVDouble("yposition"); - radSmTor = _cal->caloInfo().getDouble("radSmTor"); - xsmall = _cal->caloInfo().getDouble("radSmTor"); - xdistance = _cal->caloInfo().getDouble("xdistance"); - rInnerManifold = _cal->caloInfo().getDouble("rInnerManifold"); + phi_lbd = _cal->G4Info().get>("largeTorPhi"); + phi_sbd = _cal->G4Info().get>("smallTorPhi"); + phi_end = _cal->G4Info().get>("straightEndPhi"); + ysmall = _cal->G4Info().get>("yposition"); + radSmTor = _cal->G4Info().get("radSmTor"); + xsmall = _cal->G4Info().get("radSmTor"); + xdistance = _cal->G4Info().get("xdistance"); + rInnerManifold = _cal->G4Info().get("rInnerManifold"); } diff --git a/EventGenerator/src/CaloTBGun_module.cc b/EventGenerator/src/CaloTBGun_module.cc index e1020f0c00..04e62d124b 100644 --- a/EventGenerator/src/CaloTBGun_module.cc +++ b/EventGenerator/src/CaloTBGun_module.cc @@ -67,7 +67,7 @@ namespace mu2e { const Calorimeter& cal = *(GeomHandle()); double angleRad = angle_*M_PI/180; - double Zbuffer = cal.disk(0).crystal(400).position().z()-cal.disk(0).geomInfo().origin().z()+cal.disk(0).geomInfo().size().z()/2.0+1; + double Zbuffer = cal.disk(0).crystal(400).position().z()-cal.disk(0).diskInfo().origin().z()+cal.disk(0).diskInfo().size().z()/2.0+1; double dz = (frontVD_)? Zbuffer : 0; CLHEP::Hep3Vector pos = cal.disk(0).crystal(400).position() - CLHEP::Hep3Vector(tan(angleRad)*dz,0,dz); CLHEP::HepLorentzVector mom(energy_*sin(angleRad),0.0,energy_*cos(angleRad),energy_); diff --git a/EventGenerator/src/CosmicCRY.cc b/EventGenerator/src/CosmicCRY.cc index 583dd6a328..b584473df1 100644 --- a/EventGenerator/src/CosmicCRY.cc +++ b/EventGenerator/src/CosmicCRY.cc @@ -149,7 +149,7 @@ void CosmicCRY::generate(GenParticleCollection& genParts) { else if (_refPointChoice == "CALO") { GeomHandle calorimeter; _cosmicReferencePointInMu2e = - Hep3Vector(detsys->getOrigin().x(), _refY0, calorimeter->disk(0).geomInfo().origin().z()); + Hep3Vector(detsys->getOrigin().x(), _refY0, calorimeter->disk(0).diskInfo().origin().z()); } else if (_refPointChoice == "UNDEFINED") _cosmicReferencePointInMu2e = Hep3Vector(_refX0, _refY0, _refZ0); diff --git a/Filters/src/CaloDtsClusterFilter_module.cc b/Filters/src/CaloDtsClusterFilter_module.cc index df18cccfcf..908bad2dd4 100644 --- a/Filters/src/CaloDtsClusterFilter_module.cc +++ b/Filters/src/CaloDtsClusterFilter_module.cc @@ -20,7 +20,6 @@ // Offline #include "Offline/CalorimeterGeom/inc/Calorimeter.hh" -#include "Offline/CalorimeterGeom/inc/CaloGeomUtil.hh" #include "Offline/GeometryService/inc/GeomHandle.hh" #include "Offline/GlobalConstantsService/inc/GlobalConstantsHandle.hh" #include "Offline/GlobalConstantsService/inc/PhysicsParams.hh" @@ -141,9 +140,8 @@ namespace mu2e { CLHEP::Hep3Vector CaloDtsClusterFilter::stepPosition(const CaloShowerStep& step) const { // Get the position of the CaloShowerStep in tracker coordinates, so (x,y) are centered on (0,0) const CLHEP::Hep3Vector& pos = step.position(); - const auto& calo_geom = calorimeter_->geomUtil(); - const auto pos_in_mu2e = calo_geom.crystalToMu2e(step.volumeG4ID(), pos); - const auto pos_in_trk = calo_geom.mu2eToTracker(pos_in_mu2e); + const auto pos_in_mu2e = calorimeter_->crystalToMu2e(step.volumeG4ID(), pos); + const auto pos_in_trk = calorimeter_->mu2eToTracker(pos_in_mu2e); return pos_in_trk; } diff --git a/Filters/src/CaloPileupDtsFilter_module.cc b/Filters/src/CaloPileupDtsFilter_module.cc index 0060edfbdf..db9271a586 100644 --- a/Filters/src/CaloPileupDtsFilter_module.cc +++ b/Filters/src/CaloPileupDtsFilter_module.cc @@ -18,7 +18,6 @@ // Offline #include "Offline/CalorimeterGeom/inc/Calorimeter.hh" -#include "Offline/CalorimeterGeom/inc/CaloGeomUtil.hh" #include "Offline/DataProducts/inc/PDGCode.hh" #include "Offline/GeometryService/inc/GeomHandle.hh" #include "Offline/GlobalConstantsService/inc/GlobalConstantsHandle.hh" @@ -123,8 +122,7 @@ namespace mu2e { CLHEP::Hep3Vector CaloPileupDtsFilter::stepPosition(const CaloShowerStep& step) const { // Get the position of the CaloShowerStep in Mu2e coordinates const CLHEP::Hep3Vector& pos = step.position(); - const auto& calo_geom = calorimeter_->geomUtil(); - const auto pos_in_mu2e = calo_geom.crystalToMu2e(step.volumeG4ID(), pos); + const auto pos_in_mu2e = calorimeter_->crystalToMu2e(step.volumeG4ID(), pos); return pos_in_mu2e; } diff --git a/Filters/src/PhotonFilter_module.cc b/Filters/src/PhotonFilter_module.cc index fd04bec024..a23bbaa572 100644 --- a/Filters/src/PhotonFilter_module.cc +++ b/Filters/src/PhotonFilter_module.cc @@ -330,14 +330,14 @@ namespace mu2e { bool PhotonFilter::reportResiduals(const CaloCluster& cluster, const HelixSeedCollection& helices, const CaloClusterMC& caloClusterTruth) { - CLHEP::Hep3Vector gpos = _calorimeter->geomUtil().diskToMu2e(cluster.diskID(), cluster.cog3Vector()); + CLHEP::Hep3Vector gpos = _calorimeter->diskToMu2e(cluster.diskID(), cluster.cog3Vector()); - CLHEP::Hep3Vector tpos = _calorimeter->geomUtil().mu2eToTracker(gpos); + CLHEP::Hep3Vector tpos = _calorimeter->mu2eToTracker(gpos); - double offset = _calorimeter->caloInfo().getDouble("diskCaseZLength"); - offset += _calorimeter->caloInfo().getDouble("BPPipeZOffset"); - offset += _calorimeter->caloInfo().getDouble("BPHoleZLength"); - offset += _calorimeter->caloInfo().getDouble("FEEZLength"); + double offset = _calorimeter->G4Info().get("diskCaseZLength"); + offset += _calorimeter->G4Info().get("BPPipeZOffset"); + offset += _calorimeter->G4Info().get("BPHoleZLength"); + offset += _calorimeter->G4Info().get("FEEZLength"); offset /= 2.0; XYZVectorF clusterPos(tpos.x(), tpos.y(), tpos.z() - offset); diff --git a/GeometryService/src/DiskCalorimeterMaker.cc b/GeometryService/src/DiskCalorimeterMaker.cc index 4246efbd47..6032890336 100644 --- a/GeometryService/src/DiskCalorimeterMaker.cc +++ b/GeometryService/src/DiskCalorimeterMaker.cc @@ -53,147 +53,151 @@ namespace mu2e { verbosityLevel_ = config.getInt("calorimeter.verbosityLevel",0); - calo_->caloInfo_.set("caloDiskRadiusIn", config.getDouble("calorimeter.caloDiskRadiusIn") ); //inner radius of disk enveloppe - calo_->caloInfo_.set("caloDiskRadiusOut", config.getDouble("calorimeter.caloDiskRadiusOut") ); //outer radius of disk enveloppe - calo_->caloInfo_.set("caloFEBRadiusOut", config.getDouble("calorimeter.caloFEBRadiusOut") ); //outer radius of FEB enveloppe - calo_->caloInfo_.set("caloMotherZ0", config.getDouble("calorimeter.caloMotherZ0") ); //upstream Z of calorimeter enveloppe volume - calo_->caloInfo_.set("caloMotherZ1", config.getDouble("calorimeter.caloMotherZ1") ); //downstream Z of calorimeter enveloppe volume - calo_->caloInfo_.set("vdThickness", config.getDouble("calorimeter.vdThickness") ); //virtual detector thickness - - calo_->caloInfo_.set("diskInAlRingRIn", config.getDouble("calorimeter.diskInAlRingRIn") ); //inner radius of Al ring inside disk - calo_->caloInfo_.set("diskInAlRingZLength", config.getDouble("calorimeter.diskInAlRingZLength") ); //Z length of Al ring inside disk - calo_->caloInfo_.set("diskInCFRingRIn", config.getDouble("calorimeter.diskInCFRingRIn") ); //inner radius of carbon fiber ring - calo_->caloInfo_.set("diskInCFRingROut", config.getDouble("calorimeter.diskInCFRingROut") ); //outer radius of carbon fiber ring - calo_->caloInfo_.set("diskCaseZLength", config.getDouble("calorimeter.crystalZLength") //Z Length of crystal support structure (large Al ring) + calo_->G4Info_.set("caloDiskRadiusIn", config.getDouble("calorimeter.caloDiskRadiusIn") ); //inner radius of disk enveloppe + calo_->G4Info_.set("caloDiskRadiusOut", config.getDouble("calorimeter.caloDiskRadiusOut") ); //outer radius of disk enveloppe + calo_->G4Info_.set("caloFEBRadiusOut", config.getDouble("calorimeter.caloFEBRadiusOut") ); //outer radius of FEB enveloppe + calo_->G4Info_.set("caloMotherZ0", config.getDouble("calorimeter.caloMotherZ0") ); //upstream Z of calorimeter enveloppe volume + calo_->G4Info_.set("caloMotherZ1", config.getDouble("calorimeter.caloMotherZ1") ); //downstream Z of calorimeter enveloppe volume + calo_->G4Info_.set("caloMotherXorig", config.getDouble("mu2e.solenoidOffset") ); //x coordinate center of calo mother volume + calo_->G4Info_.set("vdThickness", config.getDouble("calorimeter.vdThickness") ); //virtual detector thickness + + calo_->G4Info_.set("diskInAlRingRIn", config.getDouble("calorimeter.diskInAlRingRIn") ); //inner radius of Al ring inside disk + calo_->G4Info_.set("diskInAlRingZLength", config.getDouble("calorimeter.diskInAlRingZLength") ); //Z length of Al ring inside disk + calo_->G4Info_.set("diskInCFRingRIn", config.getDouble("calorimeter.diskInCFRingRIn") ); //inner radius of carbon fiber ring + calo_->G4Info_.set("diskInCFRingROut", config.getDouble("calorimeter.diskInCFRingROut") ); //outer radius of carbon fiber ring + calo_->G4Info_.set("diskCaseZLength", config.getDouble("calorimeter.crystalZLength") //Z Length of crystal support structure (large Al ring) + config.getDouble("calorimeter.crystalCapZLength") ); - calo_->caloInfo_.set("diskCrystalRIn", config.getDouble("calorimeter.diskCrystalRIn") ); //inner radius of volume containing crystals - calo_->caloInfo_.set("diskCrystalROut", config.getDouble("calorimeter.diskCrystalROut") ); //outer radius of volume containing crystals - calo_->caloInfo_.set("diskCaseRingROut", config.getDouble("calorimeter.diskCaseRingROut") ); //outer radius of crystal support structure (large Al ring) - calo_->caloInfo_.set("diskOutRailZLength", config.getDouble("calorimeter.diskOutRailZLength") ); //Z length of rail outside the crystal support structure - calo_->caloInfo_.set("diskOutRailROut", config.getDouble("calorimeter.diskOutRailROut") ); //outer radius of rail outside the crystal support structure - calo_->caloInfo_.set("diskStepThickness", config.getDouble("calorimeter.diskStepThickness") ); //thickness of carbon fiber steps glued to inner hole - calo_->caloInfo_.set("diskZMotherShift", getVDouble(config,"diskZMotherShift",CaloConst::_nDisk)); //Z offset betweenfront face of disk+FEB envelope and calo mother volume Z0 - - calo_->caloInfo_.set("crystalXYLength", config.getDouble("calorimeter.crystalXYLength") ); //crysal transverse dimension - calo_->caloInfo_.set("crystalZLength", config.getDouble("calorimeter.crystalZLength") ); //crystal Z length - calo_->caloInfo_.set("crystalCapZLength", config.getDouble("calorimeter.crystalCapZLength") ); //crystal cap Z length (only in front of crystal) - calo_->caloInfo_.set("wrapperThickness", config.getDouble("calorimeter.wrapperThickness") ); //crystal wrapper thickness - calo_->caloInfo_.set("refractiveIndex", config.getDouble("calorimeter.refractiveIndex") ); //crystal refractive index - - calo_->caloInfo_.set("nSiPMPerCrystal", CaloConst::_nSiPMPerCrystal ); //number of SiPM per crystal - calo_->caloInfo_.set("readoutXLength", config.getDouble("calorimeter.readoutXLength") ); //SiPM x length - calo_->caloInfo_.set("readoutYLength", config.getDouble("calorimeter.readoutYLength") ); //SiPM Y length - calo_->caloInfo_.set("readoutZLength", config.getDouble("calorimeter.readoutZLength") ); //SiPM Z length - calo_->caloInfo_.set("shimStepsInRowId", getVInt(config,"shimStepsInRowId") ); //indices of crystal rows w.r.t center calorimeter with inner shims - calo_->caloInfo_.set("shimStepsOutRowId", getVInt(config,"shimStepsOutRowId") ); //indices of crystal rows w.r.t center calorimeter with outer shims + calo_->G4Info_.set("diskCrystalRIn", config.getDouble("calorimeter.diskCrystalRIn") ); //inner radius of volume containing crystals + calo_->G4Info_.set("diskCrystalROut", config.getDouble("calorimeter.diskCrystalROut") ); //outer radius of volume containing crystals + calo_->G4Info_.set("diskCaseRingROut", config.getDouble("calorimeter.diskCaseRingROut") ); //outer radius of crystal support structure (large Al ring) + calo_->G4Info_.set("diskOutRailZLength", config.getDouble("calorimeter.diskOutRailZLength") ); //Z length of rail outside the crystal support structure + calo_->G4Info_.set("diskOutRailROut", config.getDouble("calorimeter.diskOutRailROut") ); //outer radius of rail outside the crystal support structure + calo_->G4Info_.set("diskStepThickness", config.getDouble("calorimeter.diskStepThickness") ); //thickness of carbon fiber steps glued to inner hole + calo_->G4Info_.set("diskZMotherShift", getVDouble(config,"diskZMotherShift",CaloConst::_nDisk)); //Z offset betweenfront face of disk+FEB envelope and calo mother volume Z0 + calo_->G4Info_.set("diskCrystalFile", config.getString("calorimeter.diskCrystalFile")); //File with measured location of crystals + + calo_->G4Info_.set("crystalXYLength", config.getDouble("calorimeter.crystalXYLength") ); //crysal transverse dimension + calo_->G4Info_.set("crystalZLength", config.getDouble("calorimeter.crystalZLength") ); //crystal Z length + calo_->G4Info_.set("crystalCapZLength", config.getDouble("calorimeter.crystalCapZLength") ); //crystal cap Z length (only in front of crystal) + calo_->G4Info_.set("wrapperThickness", config.getDouble("calorimeter.wrapperThickness") ); //crystal wrapper thickness + calo_->G4Info_.set("refractiveIndex", config.getDouble("calorimeter.refractiveIndex") ); //crystal refractive index + + calo_->G4Info_.set("nSiPMPerCrystal", int(CaloConst::_nSiPMPerCrystal)); //number of SiPM per crystal + calo_->G4Info_.set("readoutXLength", config.getDouble("calorimeter.readoutXLength") ); //SiPM x length + calo_->G4Info_.set("readoutYLength", config.getDouble("calorimeter.readoutYLength") ); //SiPM Y length + calo_->G4Info_.set("readoutZLength", config.getDouble("calorimeter.readoutZLength") ); //SiPM Z length + calo_->G4Info_.set("shimStepsInRowId", getVInt(config,"shimStepsInRowId") ); //indices of crystal rows w.r.t center calorimeter with inner shims + calo_->G4Info_.set("shimStepsOutRowId", getVInt(config,"shimStepsOutRowId") ); //indices of crystal rows w.r.t center calorimeter with outer shims std::vector tempInt; for (uint16_t i=0;icaloInfo_.set("caphriCrystalId",tempInt); //LYSO crystal ids - - calo_->caloInfo_.set("hasCrates", config.getBool ("calorimeter.hasCrates") ); //include or exclude crates - calo_->caloInfo_.set("hasFrontPanel", config.getBool ("calorimeter.hasFrontPanel") ); //include or exclude front panel with calibration pipes - calo_->caloInfo_.set("hasBackPanel", config.getBool ("calorimeter.hasBackPanel") ); //include or exclude back panel with readouts and cooling pipes - calo_->caloInfo_.set("FEEXLength", config.getDouble("calorimeter.FEEXLength") ); //X length of FEE copper box - calo_->caloInfo_.set("FEEYLength", config.getDouble("calorimeter.FEEYLength") ); //Y length of FEE copper box - calo_->caloInfo_.set("FEEZLength", config.getDouble("calorimeter.FEEZLength") ); //Z length of FEE copper box - calo_->caloInfo_.set("FEEBoxThickness", config.getDouble("calorimeter.FEEBoxThickness") ); //FEE copper box thickness - calo_->caloInfo_.set("BPStripThickness", config.getDouble("calorimeter.BPStripThickness") ); //Back plate cooling strip thickness - calo_->caloInfo_.set("BPHoleXLength", config.getDouble("calorimeter.BPHoleXLength") ); //X length of hole in back plate - calo_->caloInfo_.set("BPHoleYLength", config.getDouble("calorimeter.BPHoleYLength") ); //Y length of hole in back plate - calo_->caloInfo_.set("BPHoleZLength", config.getDouble("calorimeter.BPHoleZLength") ); //Z length of hole in back plate - calo_->caloInfo_.set("BPOuterRadius", config.getDouble("calorimeter.BPOuterRadius") ); //Back plate outer radius - calo_->caloInfo_.set("BPPipeRadiusHigh", config.getDouble("calorimeter.BPPipeRadiusHigh") ); //Radius of large back plate outer pipe - calo_->caloInfo_.set("BPPipeRadiusLow", config.getDouble("calorimeter.BPPipeRadiusLow") ); //Radius of small back plate outer pipe - calo_->caloInfo_.set("BPPipeThickness", config.getDouble("calorimeter.BPPipeThickness") ); //Thickness of back plate outer pipe - calo_->caloInfo_.set("BPPipeZOffset", config.getDouble("calorimeter.BPPipeZOffset") ); //Z offset between outer pipe center and back plate - - calo_->caloInfo_.set("FPInnerRadius", config.getDouble("calorimeter.FPInnerRadius") ); //Front plate inner radius - calo_->caloInfo_.set("FPOuterRadius", config.getDouble("calorimeter.FPOuterRadius") ); //Front plate outer radius - calo_->caloInfo_.set("FPFoamZLength", config.getDouble("calorimeter.FPFoamZLength") ); //Thickness of foam inside front plate - calo_->caloInfo_.set("FPCarbonZLength", config.getDouble("calorimeter.FPCarbonZLength") ); //Z length of carbon fiber panel on front plate - calo_->caloInfo_.set("FPCoolPipeTorRadius", config.getDouble("calorimeter.FPCoolPipeTorRadius") ); //Tor radius of front plate outer cooling pipe - calo_->caloInfo_.set("FPCoolPipeRadius", config.getDouble("calorimeter.FPCoolPipeRadius") ); //Inner radius of front plate outer cooling pipe - calo_->caloInfo_.set("FPCoolPipeThickness", config.getDouble("calorimeter.FPCoolPipeThickness") ); //Thickness of front plate outer cooling pipe - calo_->caloInfo_.set("nPipes", config.getInt ("calorimeter.nPipes") ); //Number of cooling pipes inside front plate - calo_->caloInfo_.set("pipeRadius", config.getDouble("calorimeter.pipeRadius") ); //Inner radius of cooling pipes - calo_->caloInfo_.set("pipeThickness", config.getDouble("calorimeter.pipeThickness") ); //Thickness of cooling pipes - calo_->caloInfo_.set("pipeInitSeparation", config.getDouble("calorimeter.pipeInitSeparation") ); //Cooling pipe characteristics (see constructDiskCalorimeter) - calo_->caloInfo_.set("radSmTor", config.getDouble("calorimeter.radSmTor") ); //Cooling pipe characteristics (see constructDiskCalorimeter) - calo_->caloInfo_.set("xsmall", config.getDouble("calorimeter.xsmall") ); //Cooling pipe characteristics (see constructDiskCalorimeter) - calo_->caloInfo_.set("xdistance", config.getDouble("calorimeter.xdistance") ); //Cooling pipe characteristics (see constructDiskCalorimeter) - calo_->caloInfo_.set("rInnerManifold", config.getDouble("calorimeter.rInnerManifold") ); //Cooling pipe characteristics (see constructDiskCalorimeter) - calo_->caloInfo_.set("pipeTorRadius", getVDouble(config,"pipeTorRadius", calo_->caloInfo_.getInt("nPipes")) ); //Tor radius of cooling pipes and other - calo_->caloInfo_.set("largeTorPhi", getVDouble(config,"largeTorPhi", calo_->caloInfo_.getInt("nPipes")) ); //cooling pipe characteristics - calo_->caloInfo_.set("smallTorPhi", getVDouble(config,"smallTorPhi", calo_->caloInfo_.getInt("nPipes")) ); //(see constructDiskCalorimeter) - calo_->caloInfo_.set("yposition", getVDouble(config,"yposition", calo_->caloInfo_.getInt("nPipes")) ); - calo_->caloInfo_.set("straightEndPhi", getVDouble(config,"straightEndPhi",calo_->caloInfo_.getInt("nPipes")) ); - - calo_->caloInfo_.set("nCrates", config.getInt ("calorimeter.nCrates") ); //number of FEB crates per disk - calo_->caloInfo_.set("nBoards", config.getInt ("calorimeter.nBoards") ); //number of boards per crate - calo_->caloInfo_.set("FEBToDiskZOffset", config.getDouble("calorimeter.FEBToDiskZOffset") ); //Z offset between front face of FEB and front face of disk - calo_->caloInfo_.set("crateXLength", config.getDouble("calorimeter.crateXLength") ); //X length of crate - calo_->caloInfo_.set("crateYLength", config.getDouble("calorimeter.crateYLength") ); //Y length of crate - calo_->caloInfo_.set("crateZLength", config.getDouble("calorimeter.crateZLength") ); //Z length of crate - calo_->caloInfo_.set("crateFShieldThickness", config.getDouble("calorimeter.crateFShieldThickness") );//Front shield thickness - calo_->caloInfo_.set("crateBShieldThickness", config.getDouble("calorimeter.crateBShieldThickness") );//Bottom shield thickness - calo_->caloInfo_.set("crateBShieldLength", config.getDouble("calorimeter.crateBShieldLength") ); //Bottom shield length - calo_->caloInfo_.set("crateTThickness", config.getDouble("calorimeter.crateTThickness") ); //Thickness of top crate plate - calo_->caloInfo_.set("crateSThickness", config.getDouble("calorimeter.crateSThickness") ); //Thickness of side crate plate - calo_->caloInfo_.set("crateFShieldYLength", config.getDouble("calorimeter.crateFShieldYLength") ); //Y length front shield - calo_->caloInfo_.set("crateFShieldDeltaZ", config.getDouble("calorimeter.crateFShieldDeltaZ") ); //Z difference betrween front shield and crate - calo_->caloInfo_.set("radiatorThickness", config.getDouble("calorimeter.radiatorThickness") ); //radiator thickness - calo_->caloInfo_.set("radiatorZLength", config.getDouble("calorimeter.radiatorZLength") ); //radiator Z length - calo_->caloInfo_.set("activeStripThickness", config.getDouble("calorimeter.activeStripThickness") ); //active strip thickness - calo_->caloInfo_.set("passiveStripThickness", config.getDouble("calorimeter.passiveStripThickness") );//passive strip thickness - calo_->caloInfo_.set("cratePhiAngles", getVDouble(config,"cratePhiAngles",calo_->caloInfo_.getInt("nCrates")) ); //phi angle of the crates - - - // CALORIMETER ORIGIN AND FRONT FACE (FF) COORDINATES + calo_->G4Info_.set("caphriCrystalId",tempInt); //LYSO crystal ids + + calo_->G4Info_.set("hasCrates", config.getBool ("calorimeter.hasCrates") ); //include or exclude crates + calo_->G4Info_.set("hasFrontPanel", config.getBool ("calorimeter.hasFrontPanel") ); //include or exclude front panel with calibration pipes + calo_->G4Info_.set("hasBackPanel", config.getBool ("calorimeter.hasBackPanel") ); //include or exclude back panel with readouts and cooling pipes + calo_->G4Info_.set("FEEXLength", config.getDouble("calorimeter.FEEXLength") ); //X length of FEE copper box + calo_->G4Info_.set("FEEYLength", config.getDouble("calorimeter.FEEYLength") ); //Y length of FEE copper box + calo_->G4Info_.set("FEEZLength", config.getDouble("calorimeter.FEEZLength") ); //Z length of FEE copper box + calo_->G4Info_.set("FEEBoxThickness", config.getDouble("calorimeter.FEEBoxThickness") ); //FEE copper box thickness + calo_->G4Info_.set("BPStripThickness", config.getDouble("calorimeter.BPStripThickness") ); //Back plate cooling strip thickness + calo_->G4Info_.set("BPHoleXLength", config.getDouble("calorimeter.BPHoleXLength") ); //X length of hole in back plate + calo_->G4Info_.set("BPHoleYLength", config.getDouble("calorimeter.BPHoleYLength") ); //Y length of hole in back plate + calo_->G4Info_.set("BPHoleZLength", config.getDouble("calorimeter.BPHoleZLength") ); //Z length of hole in back plate + calo_->G4Info_.set("BPOuterRadius", config.getDouble("calorimeter.BPOuterRadius") ); //Back plate outer radius + calo_->G4Info_.set("BPPipeRadiusHigh", config.getDouble("calorimeter.BPPipeRadiusHigh") ); //Radius of large back plate outer pipe + calo_->G4Info_.set("BPPipeRadiusLow", config.getDouble("calorimeter.BPPipeRadiusLow") ); //Radius of small back plate outer pipe + calo_->G4Info_.set("BPPipeThickness", config.getDouble("calorimeter.BPPipeThickness") ); //Thickness of back plate outer pipe + calo_->G4Info_.set("BPPipeZOffset", config.getDouble("calorimeter.BPPipeZOffset") ); //Z offset between outer pipe center and back plate + + calo_->G4Info_.set("FPInnerRadius", config.getDouble("calorimeter.FPInnerRadius") ); //Front plate inner radius + calo_->G4Info_.set("FPOuterRadius", config.getDouble("calorimeter.FPOuterRadius") ); //Front plate outer radius + calo_->G4Info_.set("FPFoamZLength", config.getDouble("calorimeter.FPFoamZLength") ); //Thickness of foam inside front plate + calo_->G4Info_.set("FPCarbonZLength", config.getDouble("calorimeter.FPCarbonZLength") ); //Z length of carbon fiber panel on front plate + calo_->G4Info_.set("FPCoolPipeTorRadius", config.getDouble("calorimeter.FPCoolPipeTorRadius") ); //Tor radius of front plate outer cooling pipe + calo_->G4Info_.set("FPCoolPipeRadius", config.getDouble("calorimeter.FPCoolPipeRadius") ); //Inner radius of front plate outer cooling pipe + calo_->G4Info_.set("FPCoolPipeThickness", config.getDouble("calorimeter.FPCoolPipeThickness") ); //Thickness of front plate outer cooling pipe + calo_->G4Info_.set("nPipes", config.getInt ("calorimeter.nPipes") ); //Number of cooling pipes inside front plate + calo_->G4Info_.set("pipeRadius", config.getDouble("calorimeter.pipeRadius") ); //Inner radius of cooling pipes + calo_->G4Info_.set("pipeThickness", config.getDouble("calorimeter.pipeThickness") ); //Thickness of cooling pipes + calo_->G4Info_.set("pipeInitSeparation", config.getDouble("calorimeter.pipeInitSeparation") ); //Cooling pipe characteristics (see constructDiskCalorimeter) + calo_->G4Info_.set("radSmTor", config.getDouble("calorimeter.radSmTor") ); //Cooling pipe characteristics (see constructDiskCalorimeter) + calo_->G4Info_.set("xsmall", config.getDouble("calorimeter.xsmall") ); //Cooling pipe characteristics (see constructDiskCalorimeter) + calo_->G4Info_.set("xdistance", config.getDouble("calorimeter.xdistance") ); //Cooling pipe characteristics (see constructDiskCalorimeter) + calo_->G4Info_.set("rInnerManifold", config.getDouble("calorimeter.rInnerManifold") ); //Cooling pipe characteristics (see constructDiskCalorimeter) + calo_->G4Info_.set("pipeTorRadius", getVDouble(config,"pipeTorRadius", calo_->G4Info_.get("nPipes")) ); //Tor radius of cooling pipes and other + calo_->G4Info_.set("largeTorPhi", getVDouble(config,"largeTorPhi", calo_->G4Info_.get("nPipes")) ); //cooling pipe characteristics + calo_->G4Info_.set("smallTorPhi", getVDouble(config,"smallTorPhi", calo_->G4Info_.get("nPipes")) ); //(see constructDiskCalorimeter) + calo_->G4Info_.set("yposition", getVDouble(config,"yposition", calo_->G4Info_.get("nPipes")) ); + calo_->G4Info_.set("straightEndPhi", getVDouble(config,"straightEndPhi",calo_->G4Info_.get("nPipes")) ); + + calo_->G4Info_.set("nCrates", config.getInt ("calorimeter.nCrates") ); //number of FEB crates per disk + calo_->G4Info_.set("nBoards", config.getInt ("calorimeter.nBoards") ); //number of boards per crate + calo_->G4Info_.set("FEBToDiskZOffset", config.getDouble("calorimeter.FEBToDiskZOffset") ); //Z offset between front face of FEB and front face of disk + calo_->G4Info_.set("crateXLength", config.getDouble("calorimeter.crateXLength") ); //X length of crate + calo_->G4Info_.set("crateYLength", config.getDouble("calorimeter.crateYLength") ); //Y length of crate + calo_->G4Info_.set("crateZLength", config.getDouble("calorimeter.crateZLength") ); //Z length of crate + calo_->G4Info_.set("crateFShieldThickness", config.getDouble("calorimeter.crateFShieldThickness") );//Front shield thickness + calo_->G4Info_.set("crateBShieldThickness", config.getDouble("calorimeter.crateBShieldThickness") );//Bottom shield thickness + calo_->G4Info_.set("crateBShieldLength", config.getDouble("calorimeter.crateBShieldLength") ); //Bottom shield length + calo_->G4Info_.set("crateTThickness", config.getDouble("calorimeter.crateTThickness") ); //Thickness of top crate plate + calo_->G4Info_.set("crateSThickness", config.getDouble("calorimeter.crateSThickness") ); //Thickness of side crate plate + calo_->G4Info_.set("crateFShieldYLength", config.getDouble("calorimeter.crateFShieldYLength") ); //Y length front shield + calo_->G4Info_.set("crateFShieldDeltaZ", config.getDouble("calorimeter.crateFShieldDeltaZ") ); //Z difference betrween front shield and crate + calo_->G4Info_.set("radiatorThickness", config.getDouble("calorimeter.radiatorThickness") ); //radiator thickness + calo_->G4Info_.set("radiatorZLength", config.getDouble("calorimeter.radiatorZLength") ); //radiator Z length + calo_->G4Info_.set("activeStripThickness", config.getDouble("calorimeter.activeStripThickness") ); //active strip thickness + calo_->G4Info_.set("passiveStripThickness", config.getDouble("calorimeter.passiveStripThickness") );//passive strip thickness + calo_->G4Info_.set("cratePhiAngles", getVDouble(config,"cratePhiAngles",calo_->G4Info_.get("nCrates")) ); //phi angle of the crates + + + // CALORIMETER AND TRACKER ORIGIN double zTrackerCenter = config.getDouble("mu2e.detectorSystemZ0"); double xOrigin = -config.getDouble("mu2e.solenoidOffset"); double zCaloStart = config.getDouble("calorimeter.caloMotherZ0"); double zCaloEnd = config.getDouble("calorimeter.caloMotherZ1"); + calo_->G4Info_.set("caloOrigin", CLHEP::Hep3Vector(xOrigin,0,0.5*(zCaloStart+zCaloEnd)) ); //origin of the Calorimeter mother volume in Mu2e + calo_->trackerCenter_ = CLHEP::Hep3Vector(xOrigin,0,zTrackerCenter); //origin of the Tracker mother volume in Mu2e - calo_->geomUtil_.origin( CLHEP::Hep3Vector(xOrigin,0,0.5*(zCaloStart+zCaloEnd)) ); - calo_->geomUtil_.trackerCenter(CLHEP::Hep3Vector(xOrigin,0,zTrackerCenter)); - calo_->geomUtil_.crystalZLength(config.getDouble("calorimeter.crystalZLength")); - checkIt(); makeIt(); + checkIt(); } void DiskCalorimeterMaker::makeIt(void) { - double vdThickness = calo_->caloInfo_.getDouble("vdThickness"); - double caloMotherZ0 = calo_->caloInfo_.getDouble("caloMotherZ0"); - double caloMotherZ1 = calo_->caloInfo_.getDouble("caloMotherZ1"); - double crystalHalfXY = calo_->caloInfo_.getDouble("crystalXYLength")/2.0; - double wrapperHalfThick = calo_->caloInfo_.getDouble("wrapperThickness")/2.0; - double crystalHalfZLength = calo_->caloInfo_.getDouble("crystalZLength")/2.0; - double crystalCapHalfZLength = calo_->caloInfo_.getDouble("crystalCapZLength")/2.0; + double vdThickness = calo_->G4Info_.get("vdThickness"); + double caloMotherZ0 = calo_->G4Info_.get("caloMotherZ0"); + double caloMotherZ1 = calo_->G4Info_.get("caloMotherZ1"); + double crystalHalfXY = calo_->G4Info_.get("crystalXYLength")/2.0; + double wrapperHalfThick = calo_->G4Info_.get("wrapperThickness")/2.0; + double crystalHalfZLength = calo_->G4Info_.get("crystalZLength")/2.0; + double crystalCapHalfZLength = calo_->G4Info_.get("crystalCapZLength")/2.0; double crystalCellRadius = crystalHalfXY + 2.0*wrapperHalfThick; - double innerDiskRadius = calo_->caloInfo_.getDouble("diskInAlRingRIn"); - double innerCrysRadius = calo_->caloInfo_.getDouble("diskCrystalRIn"); - double outerCrysRadius = calo_->caloInfo_.getDouble("diskCrystalROut"); - double outerCaseRadius = calo_->caloInfo_.getDouble("diskCaseRingROut"); - double diskOutRailROut = calo_->caloInfo_.getDouble("diskOutRailROut"); - auto shimStepsInRowId = calo_->caloInfo_.getVInt("shimStepsInRowId"); - auto shimStepsOutRowId = calo_->caloInfo_.getVInt("shimStepsOutRowId"); - double FPCarbonThick = calo_->caloInfo_.getDouble("FPCarbonZLength"); - double FPFoamThick = calo_->caloInfo_.getDouble("FPFoamZLength"); - double FPCoolPipeRadius = calo_->caloInfo_.getDouble("FPCoolPipeRadius"); - double FPpipeRadius = calo_->caloInfo_.getDouble("pipeRadius"); + double innerDiskRadius = calo_->G4Info_.get("diskInAlRingRIn"); + double innerCrysRadius = calo_->G4Info_.get("diskCrystalRIn"); + double outerCrysRadius = calo_->G4Info_.get("diskCrystalROut"); + double outerCaseRadius = calo_->G4Info_.get("diskCaseRingROut"); + double diskOutRailROut = calo_->G4Info_.get("diskOutRailROut"); + auto shimStepsInRowId = calo_->G4Info_.get>("shimStepsInRowId"); + auto shimStepsOutRowId = calo_->G4Info_.get>("shimStepsOutRowId"); + double FPCarbonThick = calo_->G4Info_.get("FPCarbonZLength"); + double FPFoamThick = calo_->G4Info_.get("FPFoamZLength"); + double FPCoolPipeRadius = calo_->G4Info_.get("FPCoolPipeRadius"); + double FPpipeRadius = calo_->G4Info_.get("pipeRadius"); double diskCaseHalfZLength = crystalHalfZLength + crystalCapHalfZLength; - double BPHoleHalfZ = calo_->caloInfo_.getDouble("BPHoleZLength")/2.0; - double FEEBoxHalfZ = calo_->caloInfo_.getDouble("FEEZLength")/2.0; - double FEEBoxThick = calo_->caloInfo_.getDouble("FEEBoxThickness"); - double BPPipeRadiusHigh = calo_->caloInfo_.getDouble("BPPipeRadiusHigh"); - double BPPipeHalfZOffset = calo_->caloInfo_.getDouble("BPPipeZOffset")/2.0; - double crateZLength = calo_->caloInfo_.getDouble("crateZLength"); - double crateFShieldThick = calo_->caloInfo_.getDouble("crateFShieldThickness"); - double crateFShieldDeltaZ = calo_->caloInfo_.getDouble("crateFShieldDeltaZ"); - double FEBToDiskZOffset = calo_->caloInfo_.getDouble("FEBToDiskZOffset"); + double BPHoleHalfZ = calo_->G4Info_.get("BPHoleZLength")/2.0; + double FEEBoxHalfZ = calo_->G4Info_.get("FEEZLength")/2.0; + double FEEBoxThick = calo_->G4Info_.get("FEEBoxThickness"); + double BPPipeRadiusHigh = calo_->G4Info_.get("BPPipeRadiusHigh"); + double BPPipeHalfZOffset = calo_->G4Info_.get("BPPipeZOffset")/2.0; + double crateZLength = calo_->G4Info_.get("crateZLength"); + double crateFShieldThick = calo_->G4Info_.get("crateFShieldThickness"); + double crateFShieldDeltaZ = calo_->G4Info_.get("crateFShieldDeltaZ"); + double FEBToDiskZOffset = calo_->G4Info_.get("FEBToDiskZOffset"); double motherHalfZ = (caloMotherZ1 - caloMotherZ0)/2.0; + const auto& diskCrystalFile = calo_->G4Info_.get("diskCrystalFile"); + const auto& caloOrigin = calo_->G4Info_.get("caloOrigin"); + // First, calculate the total z length of the disk and the feb since we are not allowed @@ -214,9 +218,9 @@ namespace mu2e { // Then create the disks + size_t crystalOffset(0); for (int idisk=0; idiskfullCrystalList_.size(); - double separation = calo_->caloInfo_.getVDouble("diskZMotherShift").at(idisk); + double separation = calo_->G4Info_.get>("diskZMotherShift").at(idisk); double angleZ = 0; double dR1 = innerDiskRadius; @@ -226,38 +230,38 @@ namespace mu2e { CLHEP::Hep3Vector size(dR1,dR2,dZ) ; CLHEP::Hep3Vector originLocal(0, 0, -motherHalfZ + diskHalfZLength + separation + crateToDiskDeltaZ); - CLHEP::Hep3Vector frontFaceCenter = calo_->geomUtil_.origin() + originLocal + diskOriginToCrystalOrigin; + CLHEP::Hep3Vector frontFaceCenter = caloOrigin + originLocal + diskOriginToCrystalOrigin; CLHEP::Hep3Vector backFaceCenter = frontFaceCenter + CLHEP::Hep3Vector(0,0,2.0*crystalHalfZLength); CLHEP::HepRotation diskRotation = CLHEP::HepRotation::IDENTITY*CLHEP::HepRotationZ(angleZ); - auto thisDisk = std::make_shared(idisk,innerCrysRadius,outerCrysRadius, 2.0*crystalCellRadius, - 2.0*crystalHalfZLength, crystalOffset, diskOriginToCrystalOrigin); - calo_->disks_.push_back(thisDisk); - - thisDisk->geomInfo().size(size); - thisDisk->geomInfo().originLocal(originLocal); - thisDisk->geomInfo().origin(calo_->geomUtil_.origin() + originLocal); - thisDisk->geomInfo().originToCrystalOrigin(diskOriginToCrystalOrigin); - thisDisk->geomInfo().rotation(diskRotation); - thisDisk->geomInfo().frontFaceCenter(frontFaceCenter); - thisDisk->geomInfo().backFaceCenter(backFaceCenter); - thisDisk->geomInfo().FEBZOffset(crateToDiskDeltaZ); - thisDisk->geomInfo().FEBZLength(FEBZLength); - thisDisk->geomInfo().envelopeRad(dR1,dR2); - thisDisk->geomInfo().crystalDirection(CLHEP::Hep3Vector(0,0,1)); - + auto thisDisk = Disk(idisk,innerCrysRadius,outerCrysRadius, 2.0*crystalCellRadius, + 2.0*crystalHalfZLength, crystalOffset, diskOriginToCrystalOrigin, + diskCrystalFile); + + thisDisk.diskInfo().originLocal(originLocal); + thisDisk.diskInfo().setPose(caloOrigin + originLocal, diskRotation); + thisDisk.diskInfo().frontFaceCenter(frontFaceCenter); + thisDisk.diskInfo().backFaceCenter(backFaceCenter); + thisDisk.diskInfo().crystalDirection(CLHEP::Hep3Vector(0,0,1)); + thisDisk.diskInfo().size(size); + thisDisk.diskInfo().originToCrystalOrigin(diskOriginToCrystalOrigin); + thisDisk.diskInfo().FEBZOffset(crateToDiskDeltaZ); + thisDisk.diskInfo().FEBZLength(FEBZLength); + thisDisk.diskInfo().envelopeRad(dR1,dR2); + thisDisk.diskInfo().crystalZlength(2.0*crystalHalfZLength); //fill the full Crystal List / diskId (direct access for performance optimization) - for (unsigned icry=0;icrynCrystals();++icry){ - Crystal& thisCrystal = thisDisk->crystal(icry); - calo_->fullCrystalList_.push_back(&thisCrystal); + for (unsigned icry=0;icrycrystals_.push_back(&thisCrystal); //precompute the neighbors in the global frame - thisCrystal.setNeighbors(calo_->neighborsByLevel(icry + crystalOffset,1)); - thisCrystal.setNextNeighbors(calo_->neighborsByLevel(icry + crystalOffset,2)); + thisCrystal.setID(icry+crystalOffset); + thisCrystal.setNeighbors(thisDisk.neighbors(icry,1)); + thisCrystal.setNextNeighbors(thisDisk.neighbors(icry,2)); - //pre-compute the crystal position in the mu2e frame (aka global frame) - CLHEP::Hep3Vector globalPosition = thisDisk->geomInfo().origin() + thisDisk->geomInfo().inverseRotation()*(thisCrystal.localPosition()); + //pre-compute the crystal position in the mu2e (global) frame + CLHEP::Hep3Vector globalPosition = thisDisk.diskInfo().origin() + thisDisk.diskInfo().inverseRotation()*(thisCrystal.localPosition()); thisCrystal.setPosition(globalPosition); } @@ -265,40 +269,66 @@ namespace mu2e { std::vector par,stepsInX, stepsInY, stepsOutX, stepsOutY; int nrows = int(outerCaseRadius/2.0/crystalCellRadius)+1; for (int i=-nrows;iboundingBoxes(i,par); - if (par.empty()) continue; - double p0(par[0]),p1(par[1]); + auto cryList = thisDisk.idxFromRow(i); + if (cryList.empty()) continue; + + auto sortFunctor = [&thisDisk](const auto i, const auto j) + {return thisDisk.crystal(i).localPosition().x() < + thisDisk.crystal(j).localPosition().x();}; + std::sort(std::begin(cryList),std::end(cryList),sortFunctor); + + // calculate the ymin,ymax, xmin/xmax for positive and negative side + // there is a gap if the difference between the max negative and min positive > crystal size + double ymin(1e6),ymax(-1e6),xLowMin(1e6),xLowMax(-1e6),xUpMin(1e6),xUpMax(-1e6); + for (size_t ic=0;ic crystalCellRadius) { + double p4(xLowMax),p5(xUpMin); if (std::find(shimStepsInRowId.begin(),shimStepsInRowId.end(),i) != shimStepsInRowId.end()){ p4 += crystalCellRadius; p5 -= crystalCellRadius; } stepsInX.push_back(p4);stepsInX.push_back(p5); - stepsInY.push_back(par[2]);stepsInY.push_back(par[3]); + stepsInY.push_back(ymin);stepsInY.push_back(ymax); } } - calo_->caloInfo_.set("stepsInsideX",stepsInX); - calo_->caloInfo_.set("stepsInsideY",stepsInY); - calo_->caloInfo_.set("stepsOutsideX",stepsOutX); - calo_->caloInfo_.set("stepsOutsideY",stepsOutY); + calo_->G4Info_.set("stepsInsideX",stepsInX); + calo_->G4Info_.set("stepsInsideY",stepsInY); + calo_->G4Info_.set("stepsOutsideX",stepsOutX); + calo_->G4Info_.set("stepsOutsideY",stepsOutY); + + if (verbosityLevel_) std::cout<<"Constructed Disk "<id()<<": Rin="<geomInfo().innerEnvelopeR() - <<" Rout="<geomInfo().outerEnvelopeR() - <<" (X,Y,Z)="<geomInfo().origin()<<" local_(X,Y,Z)="<geomInfo().originLocal() - <<" with "<nCrystals()<<" crystals"< 1) thisDisk.print(); - if (verbosityLevel_ > 1) thisDisk->print() ; + crystalOffset += thisDisk.nCrystals(); + calo_->disks_.push_back(std::move(thisDisk)); } } @@ -308,56 +338,58 @@ namespace mu2e { { //Consistency check between disk and CaloConst content - for (const auto& diskPtr : calo_->diskPtrs()) { - if (diskPtr->nCrystals() != CaloConst::_nCrystalPerDisk) - throw cet::exception("DiskCaloGeom") << "The number of crystals ("<nCrystals() + for (size_t idisk=0;idisknDisks();++idisk) + { + const auto& thisDisk = calo_->disk(idisk); + if (thisDisk.nCrystals() != CaloConst::_nCrystalPerDisk) + throw cet::exception("DiskCaloGeom") << "The number of crystals ("<caloInfo_.getDouble("diskInAlRingRIn"); - double diskRout = calo_->caloInfo_.getDouble("diskCaseRingROut"); - double outerRingROut = calo_->caloInfo_.getDouble("diskOutRailROut"); + double diskRin = calo_->G4Info_.get("diskInAlRingRIn"); + double diskRout = calo_->G4Info_.get("diskCaseRingROut"); + double outerRingROut = calo_->G4Info_.get("diskOutRailROut"); if (diskRin > diskRout) throw cet::exception("DiskCaloGeom") << "calorimeter.diskInnerRingIn > calorimeter.diskOuterRingOut \n"; - if (outerRingROut > calo_->caloInfo_.getDouble("caloDiskRadiusOut")) + if (outerRingROut > calo_->G4Info_.get("caloDiskRadiusOut")) throw cet::exception("DiskCaloGeom") << "calorimeter outer radius larger than calorimeter mother \n"; - if (diskRin < calo_->caloInfo_.getDouble("caloDiskRadiusIn")) + if (diskRin < calo_->G4Info_.get("caloDiskRadiusIn")) throw cet::exception("DiskCaloGeom") << "calorimeter inner radius smaller than calorimeter mother \n"; //check that holes in back plate are smaller than crystal, RO smaller than holes and FEE boxes fit - if (calo_->caloInfo_.getDouble("BPHoleXLength") > calo_->caloInfo_.getDouble("crystalXYLength") || - calo_->caloInfo_.getDouble("BPHoleYLength") > calo_->caloInfo_.getDouble("crystalXYLength") ) + if (calo_->G4Info_.get("BPHoleXLength") > calo_->G4Info_.get("crystalXYLength") || + calo_->G4Info_.get("BPHoleYLength") > calo_->G4Info_.get("crystalXYLength") ) throw cet::exception("DiskCaloGeom") << "calorimeter backplate hole greater than crystal dimensions in X or Y \n"; - if (calo_->caloInfo_.getDouble("readoutXLength") > calo_->caloInfo_.getDouble("BPHoleXLength") || - calo_->caloInfo_.getDouble("readoutYLength") > calo_->caloInfo_.getDouble("BPHoleYLength")) + if (calo_->G4Info_.get("readoutXLength") > calo_->G4Info_.get("BPHoleXLength") || + calo_->G4Info_.get("readoutYLength") > calo_->G4Info_.get("BPHoleYLength")) throw cet::exception("DiskCaloGeom") << "calorimeter readout larger than hole in X or Y \n"; - if (calo_->caloInfo_.getDouble("readoutZLength") > calo_->caloInfo_.getDouble("BPHoleZLength")) + if (calo_->G4Info_.get("readoutZLength") > calo_->G4Info_.get("BPHoleZLength")) throw cet::exception("DiskCaloGeom") << "calorimeter readout too thick to fit in hole \n"; - if (calo_->caloInfo_.getDouble("FEEXLength") > calo_->caloInfo_.getDouble("BPHoleXLength")/calo_->caloInfo_.getInt("nSiPMPerCrystal")) + if (calo_->G4Info_.get("FEEXLength") > calo_->G4Info_.get("BPHoleXLength")/calo_->G4Info_.get("nSiPMPerCrystal")) throw cet::exception("DiskCaloGeom") << "calorimeter FEE box does not fit in X direction \n"; - if (calo_->caloInfo_.getDouble("FEEYLength") > calo_->caloInfo_.getDouble("crystalXYLength")-2*calo_->caloInfo_.getDouble("FEEBoxThickness")) + if (calo_->G4Info_.get("FEEYLength") > calo_->G4Info_.get("crystalXYLength")-2*calo_->G4Info_.get("FEEBoxThickness")) throw cet::exception("DiskCaloGeom") << "calorimeter FEE box does not fit in Y direction \n"; //Check pipes - for (int i=0;icaloInfo().getInt("nPipes");++i){ - double pipeTorRadius = calo_->caloInfo_.getVDouble("pipeTorRadius").at(i); - double pipeRadius = calo_->caloInfo_.getDouble("pipeRadius"); - double radiusIn = calo_->caloInfo_.getDouble("FPInnerRadius"); - double radiusOut = calo_->caloInfo_.getDouble("FPOuterRadius"); - double foamZLength = calo_->caloInfo_.getDouble("FPFoamZLength"); - double carbonThick = calo_->caloInfo_.getDouble("FPCarbonZLength"); - double foamThick = calo_->caloInfo_.getDouble("FPFoamZLength"); - double coolPipeTorRadius = calo_->caloInfo_.getDouble("FPCoolPipeTorRadius"); - double coolPipeRadius = calo_->caloInfo_.getDouble("FPCoolPipeRadius"); + for (int i=0;iG4Info().get("nPipes");++i){ + double pipeTorRadius = calo_->G4Info_.get>("pipeTorRadius").at(i); + double pipeRadius = calo_->G4Info_.get("pipeRadius"); + double radiusIn = calo_->G4Info_.get("FPInnerRadius"); + double radiusOut = calo_->G4Info_.get("FPOuterRadius"); + double foamZLength = calo_->G4Info_.get("FPFoamZLength"); + double carbonThick = calo_->G4Info_.get("FPCarbonZLength"); + double foamThick = calo_->G4Info_.get("FPFoamZLength"); + double coolPipeTorRadius = calo_->G4Info_.get("FPCoolPipeTorRadius"); + double coolPipeRadius = calo_->G4Info_.get("FPCoolPipeRadius"); if ( pipeTorRadius - pipeRadius < radiusIn) throw cet::exception("DiskCaloGeom") << "element "<caloInfo_.getInt("nBoards"); - double radiatorDY = calo_->caloInfo_.getDouble("radiatorThickness")/2.0; - double activeStripDY = calo_->caloInfo_.getDouble("activeStripThickness")/2.0; - double passiveStripDY = calo_->caloInfo_.getDouble("passiveStripThickness")/2.0; - double crateXLength = calo_->caloInfo_.getDouble("crateXLength"); - double crateYLength = calo_->caloInfo_.getDouble("crateYLength"); + int nBoards = calo_->G4Info_.get("nBoards"); + double radiatorDY = calo_->G4Info_.get("radiatorThickness")/2.0; + double activeStripDY = calo_->G4Info_.get("activeStripThickness")/2.0; + double passiveStripDY = calo_->G4Info_.get("passiveStripThickness")/2.0; + double crateXLength = calo_->G4Info_.get("crateXLength"); + double crateYLength = calo_->G4Info_.get("crateYLength"); if ( nBoards*(radiatorDY+activeStripDY+passiveStripDY) > crateYLength) throw cet::exception("DiskCaloGeom") << "calorimeter FEB boards too thick\n"; - if (calo_->caloInfo_.getDouble("crateFShieldYLength") > crateYLength) + if (calo_->G4Info_.get("crateFShieldYLength") > crateYLength) throw cet::exception("DiskCaloGeom") << "calorimeter FEB front shile too long in Y direction\n"; - if (calo_->caloInfo_.getDouble("crateSThickness") > crateXLength) + if (calo_->G4Info_.get("crateSThickness") > crateXLength) throw cet::exception("DiskCaloGeom") << "calorimeter FEB crate side too thick\n"; - if (calo_->caloInfo_.getDouble("crateTThickness") > crateYLength ) + if (calo_->G4Info_.get("crateTThickness") > crateYLength ) throw cet::exception("DiskCaloGeom") << "calorimeter FEB crate top too thick\n"; } diff --git a/GeometryService/src/VirtualDetectorMaker.cc b/GeometryService/src/VirtualDetectorMaker.cc index d036ea21ef..1a6bf0999f 100644 --- a/GeometryService/src/VirtualDetectorMaker.cc +++ b/GeometryService/src/VirtualDetectorMaker.cc @@ -435,20 +435,21 @@ namespace mu2e { int vdIdFEBEdge = VirtualDetectorId::EMC_FEB_0_EdgeIn; int vdIdFEBSurf = VirtualDetectorId::EMC_FEB_0_SurfIn; - CLHEP::Hep3Vector parentInMu2e = cg->geomUtil().origin(); + CLHEP::Hep3Vector parentInMu2e = cg->G4Info().get("caloOrigin"); - for (const auto& diskPtr : cg->diskPtrs()) + for (size_t idisk=0;idisknDisks();++idisk) { - const CLHEP::Hep3Vector& sizeDisk = diskPtr->geomInfo().size(); - CLHEP::Hep3Vector posDiskLocal = diskPtr->geomInfo().origin() - cg->geomUtil().origin(); - CLHEP::Hep3Vector posCrateLocal = posDiskLocal + CLHEP::Hep3Vector(0.0,0.0,diskPtr->geomInfo().FEBZOffset()); + const auto& thisDisk = cg->disk(idisk); + const CLHEP::Hep3Vector& sizeDisk = thisDisk.diskInfo().size(); + CLHEP::Hep3Vector posDiskLocal = thisDisk.diskInfo().origin() - parentInMu2e; + CLHEP::Hep3Vector posCrateLocal = posDiskLocal + CLHEP::Hep3Vector(0.0,0.0,thisDisk.diskInfo().FEBZOffset()); CLHEP::Hep3Vector posFrontDisk = posDiskLocal - CLHEP::Hep3Vector (0,0,sizeDisk.z()/2.0-vdHL); CLHEP::Hep3Vector posBackDisk = posDiskLocal + CLHEP::Hep3Vector (0,0,sizeDisk.z()/2.0-vdHL); CLHEP::Hep3Vector posInnerDisk = posDiskLocal; - CLHEP::Hep3Vector posFrontFEB = posFrontDisk - CLHEP::Hep3Vector(0.0,0.0,diskPtr->geomInfo().FEBZOffset()); - CLHEP::Hep3Vector posBackFEB = posFrontFEB + CLHEP::Hep3Vector (0,0,diskPtr->geomInfo().FEBZLength()); + CLHEP::Hep3Vector posFrontFEB = posFrontDisk - CLHEP::Hep3Vector(0.0,0.0,thisDisk.diskInfo().FEBZOffset()); + CLHEP::Hep3Vector posBackFEB = posFrontFEB + CLHEP::Hep3Vector (0,0,thisDisk.diskInfo().FEBZLength()); CLHEP::Hep3Vector posInnerFEB = (posBackFEB+posFrontFEB)/2.0; vd->addVirtualDetector(vdIdDiskSurf, parentInMu2e,0,posFrontDisk); diff --git a/Mu2eG4/geom/calorimeter_CsI_v2.txt b/Mu2eG4/geom/calorimeter_CsI_v2.txt index 94d925132e..27665781a6 100644 --- a/Mu2eG4/geom/calorimeter_CsI_v2.txt +++ b/Mu2eG4/geom/calorimeter_CsI_v2.txt @@ -20,6 +20,7 @@ double calorimeter.diskOutRailZLength = 30.; double calorimeter.diskOutRailROut = 719.; double calorimeter.diskStepThickness = 3.; vector calorimeter.diskZMotherShift = {0.0,700.0}; +string calorimeter.diskCrystalFile = "Offline/CalorimeterGeom/data/crystalPos.txt"; double calorimeter.crystalXYLength = 34.0; double calorimeter.crystalZLength = 200.0; diff --git a/Mu2eG4/src/CaloReadoutCardSD.cc b/Mu2eG4/src/CaloReadoutCardSD.cc index 862910ce00..4e041d89ca 100644 --- a/Mu2eG4/src/CaloReadoutCardSD.cc +++ b/Mu2eG4/src/CaloReadoutCardSD.cc @@ -28,7 +28,7 @@ namespace mu2e { Mu2eG4SensitiveDetector(name,config),_nro(0) { GeomHandle cg; - _nro = cg->caloInfo().getInt("nSiPMPerCrystal"); + _nro = cg->G4Info().get("nSiPMPerCrystal"); } G4bool CaloReadoutCardSD::ProcessHits(G4Step* aStep,G4TouchableHistory*) diff --git a/Mu2eG4/src/CaloReadoutSD.cc b/Mu2eG4/src/CaloReadoutSD.cc index fb0ea211b7..f4a2000312 100644 --- a/Mu2eG4/src/CaloReadoutSD.cc +++ b/Mu2eG4/src/CaloReadoutSD.cc @@ -29,7 +29,7 @@ namespace mu2e { Mu2eG4SensitiveDetector(name,config),_nro(0) { GeomHandle cg; - _nro = cg->caloInfo().getInt("nSiPMPerCrystal"); + _nro = cg->G4Info().get("nSiPMPerCrystal"); } G4bool CaloReadoutSD::ProcessHits(G4Step* aStep,G4TouchableHistory*) diff --git a/Mu2eG4/src/constructDiskCalorimeter.cc b/Mu2eG4/src/constructDiskCalorimeter.cc index e644d5e577..13a1485332 100644 --- a/Mu2eG4/src/constructDiskCalorimeter.cc +++ b/Mu2eG4/src/constructDiskCalorimeter.cc @@ -67,16 +67,16 @@ namespace mu2e { const int verbosity = config.getInt("calorimeter.verbosityLevel",1); const unsigned nDisks = cal.nDisks(); - const double caloDiskRadiusIn = cal.caloInfo().getDouble("caloDiskRadiusIn"); - const double caloDiskRadiusOut = cal.caloInfo().getDouble("caloDiskRadiusOut"); - const double caloFEBRadiusOut = cal.caloInfo().getDouble("caloFEBRadiusOut"); - const double mother_z0 = cal.caloInfo().getDouble("caloMotherZ0"); - const double mother_z1 = cal.caloInfo().getDouble("caloMotherZ1"); + const double caloDiskRadiusIn = cal.G4Info().get("caloDiskRadiusIn"); + const double caloDiskRadiusOut = cal.G4Info().get("caloDiskRadiusOut"); + const double caloFEBRadiusOut = cal.G4Info().get("caloFEBRadiusOut"); + const double mother_z0 = cal.G4Info().get("caloMotherZ0"); + const double mother_z1 = cal.G4Info().get("caloMotherZ1"); const double mother_zlength = mother_z1-mother_z0; const double mother_zCenter = (mother_z1+mother_z0)/2.0; - const double FEBOffset = cal.caloInfo().getDouble("FEBToDiskZOffset"); + const double FEBOffset = cal.G4Info().get("FEBToDiskZOffset"); const auto FEBPhiMinMax = calcFEBPhiRange(cal); - const bool hasCrates = cal.caloInfo().getBool("hasCrates"); + const bool hasCrates = cal.G4Info().get("hasCrates"); const bool hasCable = ds.hasCableRunCal() && hasCrates; @@ -115,10 +115,10 @@ namespace mu2e { // Disk and FEB volumes std::vector calorimeterDisk(nDisks), calorimeterFEB(nDisks); for (size_t idisk=0;idiskdoSurfaceCheck("calorimeterEnvelope"); const int verbosity = config.getInt("calorimeter.verbosityLevel",1); - const bool hasFrontPanel = cal.caloInfo().getBool("hasFrontPanel"); - const bool hasBackPanel = cal.caloInfo().getBool("hasBackPanel"); - const double vdThickness = cal.caloInfo().getDouble("vdThickness"); - const double R0disk = cal.caloInfo().getDouble("caloDiskRadiusIn"); - const double R1disk = cal.caloInfo().getDouble("caloDiskRadiusOut"); + const bool hasFrontPanel = cal.G4Info().get("hasFrontPanel"); + const bool hasBackPanel = cal.G4Info().get("hasBackPanel"); + const double vdThickness = cal.G4Info().get("vdThickness"); + const double R0disk = cal.G4Info().get("caloDiskRadiusIn"); + const double R1disk = cal.G4Info().get("caloDiskRadiusOut"); G4Material* vacuumMaterial = materialFinder.get("calorimeter.vacuumMaterial"); @@ -241,7 +241,7 @@ namespace mu2e { auto volumePtr = findCaloSolid(fullDisk.logical,"CaloCrystalCsI_"+std::to_string(idisk), volumeNodes); volumeNodes.push_back(fullDisk.logical); - //navigate the nodes to calculate the total translation, and add half crystal size to match the originToCrystalOrigin value in geomInfo + //navigate the nodes to calculate the total translation, and add half crystal size to match the originToCrystalOrigin value in diskInfo double zTranslation(0); for (size_t i=1;iGetNoDaughters();++j) { @@ -253,18 +253,18 @@ namespace mu2e { } if (volumePtr) zTranslation -= dynamic_cast(volumePtr->GetSolid())->GetZHalfLength(); - double delta1 = 2*zHalftot-cal.disk(idisk).geomInfo().size().z(); - double delta2 = zTranslation-cal.disk(idisk).geomInfo().originToCrystalOrigin().z(); + double delta1 = 2*zHalftot-cal.disk(idisk).diskInfo().size().z(); + double delta2 = zTranslation-cal.disk(idisk).diskInfo().originToCrystalOrigin().z(); if (std::abs(delta1) > 1e-3) G4cout << __func__ <<"PANIC..... geometry description in Geant4 and DiskMaker do NOT match - disk size: " - << 2*zHalftot<<" vs "< 1e-3) G4cout << __func__ <<"PANIC..... geometry description in Geant4 and DiskMaker do NOT match - originToCrystalOrigin: " - << zTranslation<<" vs "<doSurfaceCheck ("calorimeterPipe"); const int verbosity = config.getInt ("calorimeter.verbosityLevel",1); - const double FPInnerRadius = cal.caloInfo().getDouble("FPInnerRadius"); - const double FPOuterRadius = cal.caloInfo().getDouble("FPOuterRadius"); - const double FPCarbonDZ = cal.caloInfo().getDouble("FPCarbonZLength")/2.0; - const double FPFoamDZ = cal.caloInfo().getDouble("FPFoamZLength")/2.0; - const double FPCoolPipeTorRadius = cal.caloInfo().getDouble("FPCoolPipeTorRadius"); - const double FPCoolPipeRadius = cal.caloInfo().getDouble("FPCoolPipeRadius"); - const double FPCoolPipeThickness = cal.caloInfo().getDouble("FPCoolPipeThickness"); + const double FPInnerRadius = cal.G4Info().get("FPInnerRadius"); + const double FPOuterRadius = cal.G4Info().get("FPOuterRadius"); + const double FPCarbonDZ = cal.G4Info().get("FPCarbonZLength")/2.0; + const double FPFoamDZ = cal.G4Info().get("FPFoamZLength")/2.0; + const double FPCoolPipeTorRadius = cal.G4Info().get("FPCoolPipeTorRadius"); + const double FPCoolPipeRadius = cal.G4Info().get("FPCoolPipeRadius"); + const double FPCoolPipeThickness = cal.G4Info().get("FPCoolPipeThickness"); const double FPCoolPipeRadiusIn = FPCoolPipeRadius-FPCoolPipeThickness; - const int nPipes = cal.caloInfo().getInt ("nPipes"); - const double pipeRadius = cal.caloInfo().getDouble ("pipeRadius"); - const double pipeThickness = cal.caloInfo().getDouble ("pipeThickness"); - const double pipeInitSeparation = cal.caloInfo().getDouble ("pipeInitSeparation"); - const std::vector pipeTorRadius = cal.caloInfo().getVDouble("pipeTorRadius"); - const std::vector largeTorPhi = cal.caloInfo().getVDouble("largeTorPhi"); - const std::vector smallTorPhi = cal.caloInfo().getVDouble("smallTorPhi"); - const std::vector yposition = cal.caloInfo().getVDouble("yposition"); - const std::vector straightEndPhi = cal.caloInfo().getVDouble("straightEndPhi"); - const double radSmTor = cal.caloInfo().getDouble ("radSmTor"); - const double xsmall = cal.caloInfo().getDouble ("xsmall"); - const double xdistance = cal.caloInfo().getDouble ("xdistance"); + const int nPipes = cal.G4Info().get ("nPipes"); + const double pipeRadius = cal.G4Info().get ("pipeRadius"); + const double pipeThickness = cal.G4Info().get ("pipeThickness"); + const double pipeInitSeparation = cal.G4Info().get ("pipeInitSeparation"); + const std::vector pipeTorRadius = cal.G4Info().get>("pipeTorRadius"); + const std::vector largeTorPhi = cal.G4Info().get>("largeTorPhi"); + const std::vector smallTorPhi = cal.G4Info().get>("smallTorPhi"); + const std::vector yposition = cal.G4Info().get>("yposition"); + const std::vector straightEndPhi = cal.G4Info().get>("straightEndPhi"); + const double radSmTor = cal.G4Info().get ("radSmTor"); + const double xsmall = cal.G4Info().get ("xsmall"); + const double xdistance = cal.G4Info().get ("xdistance"); const double frontPanelHalfThick = (2.0*FPCarbonDZ+2.0*FPFoamDZ-pipeRadius+FPCoolPipeRadius)/2.0; const double ZposCarbon2 = frontPanelHalfThick-FPCarbonDZ; @@ -475,33 +475,33 @@ namespace mu2e { G4Material* outerRingMaterial = materialFinder.get("calorimeter.outerRingMaterial"); G4Material* coolPipeMaterial = materialFinder.get("calorimeter.coolPipeMaterial"); - const double vdThickness = cal.caloInfo().getDouble("vdThickness"); - const double crystalDXY = cal.caloInfo().getDouble("crystalXYLength")/2.0; - const double crystalDZ = cal.caloInfo().getDouble("crystalZLength")/2.0; - const double crystalCapDZ = cal.caloInfo().getDouble("crystalCapZLength")/2.0; - const double wrapperHalfThick = cal.caloInfo().getDouble("wrapperThickness")/2.0; + const double vdThickness = cal.G4Info().get("vdThickness"); + const double crystalDXY = cal.G4Info().get("crystalXYLength")/2.0; + const double crystalDZ = cal.G4Info().get("crystalZLength")/2.0; + const double crystalCapDZ = cal.G4Info().get("crystalCapZLength")/2.0; + const double wrapperHalfThick = cal.G4Info().get("wrapperThickness")/2.0; const double wrapperDXY = crystalDXY + 2*wrapperHalfThick; const double wrapperDZ = crystalDZ+crystalCapDZ; - const std::vector caphriCystalId = cal.caloInfo().getVInt("caphriCrystalId"); + const std::vector caphriCystalId = cal.G4Info().get>("caphriCrystalId"); const double diskCaseDZ = wrapperDZ; - const double diskInAlRingRIn = cal.caloInfo().getDouble("diskInAlRingRIn"); - const double diskInAlRingDZ = cal.caloInfo().getDouble("diskInAlRingZLength")/2.0; - const double diskInCFRingRIn = cal.caloInfo().getDouble("diskInCFRingRIn"); - const double diskInCFRingROut = cal.caloInfo().getDouble("diskInCFRingROut"); - const double diskCaseRingROut = cal.caloInfo().getDouble("diskCaseRingROut"); - const double diskOutRailROut = cal.caloInfo().getDouble("diskOutRailROut") - vdThickness; - const double diskOutRailDZ = cal.caloInfo().getDouble("diskOutRailZLength")/2.0; - - const double FPCoolPipeRadius = cal.caloInfo().getDouble("FPCoolPipeRadius"); - const double FPCoolPipeThickness = cal.caloInfo().getDouble("FPCoolPipeThickness"); + const double diskInAlRingRIn = cal.G4Info().get("diskInAlRingRIn"); + const double diskInAlRingDZ = cal.G4Info().get("diskInAlRingZLength")/2.0; + const double diskInCFRingRIn = cal.G4Info().get("diskInCFRingRIn"); + const double diskInCFRingROut = cal.G4Info().get("diskInCFRingROut"); + const double diskCaseRingROut = cal.G4Info().get("diskCaseRingROut"); + const double diskOutRailROut = cal.G4Info().get("diskOutRailROut") - vdThickness; + const double diskOutRailDZ = cal.G4Info().get("diskOutRailZLength")/2.0; + + const double FPCoolPipeRadius = cal.G4Info().get("FPCoolPipeRadius"); + const double FPCoolPipeThickness = cal.G4Info().get("FPCoolPipeThickness"); const double coolPipeZpos = (diskCaseDZ - 2*FPCoolPipeRadius - 2.0*diskOutRailDZ)/2.0 + FPCoolPipeRadius; - const std::vector stepsInX = cal.caloInfo().getVDouble("stepsInsideX"); - const std::vector stepsInY = cal.caloInfo().getVDouble("stepsInsideY"); - const std::vector stepsOutX = cal.caloInfo().getVDouble("stepsOutsideX"); - const std::vector stepsOutY = cal.caloInfo().getVDouble("stepsOutsideY"); - const double diskStepThickness = cal.caloInfo().getDouble ("diskStepThickness"); + const std::vector stepsInX = cal.G4Info().get>("stepsInsideX"); + const std::vector stepsInY = cal.G4Info().get>("stepsInsideY"); + const std::vector stepsOutX = cal.G4Info().get>("stepsOutsideX"); + const std::vector stepsOutY = cal.G4Info().get>("stepsOutsideY"); + const double diskStepThickness = cal.G4Info().get ("diskStepThickness"); //------------------------------------------------------------ @@ -694,36 +694,36 @@ namespace mu2e { G4Material* pipeMaterial = materialFinder.get("calorimeter.coolPipeMaterial"); G4Material* stripMaterial = materialFinder.get("calorimeter.BPStripMaterial"); - const double BPInnerRadius = cal.caloInfo().getDouble("FPInnerRadius"); //same as front plate - const double BPOuterRadius = cal.caloInfo().getDouble("BPOuterRadius"); - const double diskCrystalRIn = cal.caloInfo().getDouble("diskCrystalRIn"); - const double diskCrystalROut = cal.caloInfo().getDouble("diskCrystalROut"); - - const double crystalDXY = cal.caloInfo().getDouble("crystalXYLength")/2.0; - const double wrapperDXY = crystalDXY + cal.caloInfo().getDouble("wrapperThickness"); - const double RODX = cal.caloInfo().getDouble("readoutXLength")/2.0; - const double RODY = cal.caloInfo().getDouble("readoutYLength")/2.0; - const double RODZ = cal.caloInfo().getDouble("readoutZLength")/2.0; - const double holeDX = cal.caloInfo().getDouble("BPHoleXLength")/2.0; - const double holeDY = cal.caloInfo().getDouble("BPHoleYLength")/2.0; - const double holeDZ = cal.caloInfo().getDouble("BPHoleZLength")/2.0; + const double BPInnerRadius = cal.G4Info().get("FPInnerRadius"); //same as front plate + const double BPOuterRadius = cal.G4Info().get("BPOuterRadius"); + const double diskCrystalRIn = cal.G4Info().get("diskCrystalRIn"); + const double diskCrystalROut = cal.G4Info().get("diskCrystalROut"); + + const double crystalDXY = cal.G4Info().get("crystalXYLength")/2.0; + const double wrapperDXY = crystalDXY + cal.G4Info().get("wrapperThickness"); + const double RODX = cal.G4Info().get("readoutXLength")/2.0; + const double RODY = cal.G4Info().get("readoutYLength")/2.0; + const double RODZ = cal.G4Info().get("readoutZLength")/2.0; + const double holeDX = cal.G4Info().get("BPHoleXLength")/2.0; + const double holeDY = cal.G4Info().get("BPHoleYLength")/2.0; + const double holeDZ = cal.G4Info().get("BPHoleZLength")/2.0; const double stripDY = wrapperDXY-holeDY-1.0; - const double stripDZ = cal.caloInfo().getDouble("BPStripThickness")/2.0; - const double FEEDX = cal.caloInfo().getDouble("FEEXLength")/2.0; - const double FEEDY = cal.caloInfo().getDouble("FEEYLength")/2.0; - const double FEEDZ = cal.caloInfo().getDouble("FEEZLength")/2.0; + const double stripDZ = cal.G4Info().get("BPStripThickness")/2.0; + const double FEEDX = cal.G4Info().get("FEEXLength")/2.0; + const double FEEDY = cal.G4Info().get("FEEYLength")/2.0; + const double FEEDZ = cal.G4Info().get("FEEZLength")/2.0; - const double FEEBoxThickness = cal.caloInfo().getDouble("FEEBoxThickness"); + const double FEEBoxThickness = cal.G4Info().get("FEEBoxThickness"); const double FEEBoxDX = holeDX + 2*FEEBoxThickness; const double FEEBoxDY = FEEDY + 2*FEEBoxThickness; const double FEEBoxDZ = FEEDZ + 2*FEEBoxThickness; - const double BPPipeRadiusHigh = cal.caloInfo().getDouble("BPPipeRadiusHigh"); - const double BPPipeRadiusLow = cal.caloInfo().getDouble("BPPipeRadiusLow"); - const double BPPipeThickness = cal.caloInfo().getDouble("BPPipeThickness"); + const double BPPipeRadiusHigh = cal.G4Info().get("BPPipeRadiusHigh"); + const double BPPipeRadiusLow = cal.G4Info().get("BPPipeRadiusLow"); + const double BPPipeThickness = cal.G4Info().get("BPPipeThickness"); const double BPPipeTorRadiusHigh = BPOuterRadius - BPPipeRadiusHigh; const double BPPipeTorRadiusLow = BPOuterRadius - 3.0*BPPipeRadiusHigh; - const double BPPipeDZOffset = cal.caloInfo().getDouble("BPPipeZOffset")/2.0; + const double BPPipeDZOffset = cal.G4Info().get("BPPipeZOffset")/2.0; const double BPFEEDZ = FEEBoxDZ + BPPipeDZOffset + BPPipeRadiusHigh; @@ -770,7 +770,7 @@ namespace mu2e { for(unsigned ic=0; ic doSurfaceCheck("calorimeterCrate"); const int verbosity = config.getInt("calorimeter.verbosityLevel",1); - const int nCrates = cal.caloInfo().getInt("nCrates"); - const double vdThickness = cal.caloInfo().getDouble("vdThickness"); - const double caloDiskRadiusOut = cal.caloInfo().getDouble("caloDiskRadiusOut"); - const auto cratePhiAngle = cal.caloInfo().getVDouble("cratePhiAngles"); + const int nCrates = cal.G4Info().get("nCrates"); + const double vdThickness = cal.G4Info().get("vdThickness"); + const double caloDiskRadiusOut = cal.G4Info().get("caloDiskRadiusOut"); + const auto cratePhiAngle = cal.G4Info().get>("cratePhiAngles"); // get me a lil' crate VolumeInfo crateBox = caloBuildCrate(config,idisk); @@ -986,14 +986,14 @@ namespace mu2e { } } - double delta = dynamic_cast(fullFEB.solid)->GetZHalfLength()*2 - cal.disk(idisk).geomInfo().FEBZLength(); + double delta = dynamic_cast(fullFEB.solid)->GetZHalfLength()*2 - cal.disk(idisk).diskInfo().FEBZLength(); if (std::abs(delta) > 1e-3) G4cout << __func__ <<"PANIC..... geometry description in Geant4 and DiskMaker do NOT match - FEB size: " <(fullFEB.solid)->GetZHalfLength()*2<<" vs " - <(fullFEB.solid)->GetZHalfLength()*2 - <<" / "< ("nBoards"); + const double crateDX = cal.G4Info().get("crateXLength")/2.0; + const double crateDY = cal.G4Info().get("crateYLength")/2.0; + const double crateDZ = cal.G4Info().get("crateZLength")/2.0; + const double crateFShieldDisp = cal.G4Info().get("crateFShieldDeltaZ"); + const double crateFShieldThick = cal.G4Info().get("crateFShieldThickness"); + const double crateBottomThick = cal.G4Info().get("crateBShieldThickness"); + const double crateBottomLength = cal.G4Info().get("crateBShieldLength"); + const double crateTopThick = cal.G4Info().get("crateTThickness"); + const double crateSideThick = cal.G4Info().get("crateSThickness"); + const double crateFShieldDY = cal.G4Info().get("crateFShieldYLength")/2.0; + + const double radiatorDY = cal.G4Info().get("radiatorThickness")/2.0; + const double radiatorDZ = cal.G4Info().get("radiatorZLength")/2.0; + const double activeStripDY = cal.G4Info().get("activeStripThickness")/2.0; + const double passiveStripDY = cal.G4Info().get("passiveStripThickness")/2.0; const double crateFullDZ = crateDZ+crateFShieldDisp/2.0+crateFShieldThick/2.0; const double crateFullDY = crateDY+crateBottomThick/2.0; const double crateBoxInDY = crateDY-crateTopThick/2.0; @@ -1151,8 +1151,8 @@ namespace mu2e { G4Material* cableMaterial = findMaterialOrThrow(ds.calCableRunMaterial()); G4Material* cableCoreMaterial = findMaterialOrThrow(ds.materialCableRunCalCore()); - const double mother_z0 = cal.caloInfo().getDouble("caloMotherZ0"); - const double mother_z1 = cal.caloInfo().getDouble("caloMotherZ1"); + const double mother_z0 = cal.G4Info().get("caloMotherZ0"); + const double mother_z1 = cal.G4Info().get("caloMotherZ1"); const double mother_zlength = mother_z1-mother_z0; const double crRin = ds.upRInCableRunCal(); const double crRout = ds.upROutCableRunCal(); @@ -1162,7 +1162,7 @@ namespace mu2e { //length of the cable run = distance between disk for first gaps, ditance betweenm end of feb and mother volume for last gap G4Tubs* FEB = dynamic_cast(FEBvol.solid); - G4double distFEBtoNext = (idisk+1 < cal.nDisks()) ? cal.disk(idisk+1).geomInfo().origin().z() - cal.disk(idisk).geomInfo().origin().z() + G4double distFEBtoNext = (idisk+1 < cal.nDisks()) ? cal.disk(idisk+1).diskInfo().origin().z() - cal.disk(idisk).diskInfo().origin().z() : mother_zlength/2.0 - FEBvol.physical->GetTranslation().z()-FEB->GetZHalfLength(); G4double cableHalfLength = (idisk+1 < cal.nDisks()) ? 0.5*distFEBtoNext - FEB->GetZHalfLength() : 0.5*distFEBtoNext; @@ -1255,9 +1255,9 @@ namespace mu2e { // utility to make list of vertexes for extruded polygon solid std::vector calcFEBPhiRange(const DiskCalorimeter& cal) { - G4double crateXLength = cal.caloInfo().getDouble("crateXLength"); - G4double crateRin = cal.caloInfo().getDouble("caloDiskRadiusOut"); - auto cratePhis = cal.caloInfo().getVDouble("cratePhiAngles"); + G4double crateXLength = cal.G4Info().get("crateXLength"); + G4double crateRin = cal.G4Info().get("caloDiskRadiusOut"); + auto cratePhis = cal.G4Info().get>("cratePhiAngles"); for (auto& val : cratePhis) val *= CLHEP::degree; G4double cratePhiSpan = atan(crateXLength/2.0/crateRin); diff --git a/Mu2eKinKal/inc/KKFit.hh b/Mu2eKinKal/inc/KKFit.hh index 3b4dee0c0e..1112912c42 100644 --- a/Mu2eKinKal/inc/KKFit.hh +++ b/Mu2eKinKal/inc/KKFit.hh @@ -345,10 +345,10 @@ namespace mu2e { template std::shared_ptr KKFit::caloAxis(CaloCluster const& cluster, Calorimeter const& calo) const { // move cluster COG into the tracker frame. COG is at the front face of the disk - CLHEP::Hep3Vector cog = calo.geomUtil().mu2eToTracker(calo.geomUtil().diskFFToMu2e( cluster.diskID(), cluster.cog3Vector())); + CLHEP::Hep3Vector cog = calo.mu2eToTracker(calo.diskFFToMu2e( cluster.diskID(), cluster.cog3Vector())); // project this along the crystal axis to the SIPM, which is at the back. This is the point the time measurement corresponds to VEC3 ffcog(cog); - double lcrystal = calo.caloInfo().getDouble("crystalZLength"); // text-keyed lookup is very inefficient FIXME! + double lcrystal = calo.G4Info().get("crystalZLength"); // text-keyed lookup is very inefficient FIXME! VEC3 crystalF2B = VEC3(0.0,0.0,lcrystal); // this should come directly from the calogeometry, TODO VEC3 sipmcog = ffcog + crystalF2B; // create the SensorLine trajectory from this information: signal goes towards the sipm @@ -580,7 +580,7 @@ namespace mu2e { template void KKFit::addCaloHit(Calorimeter const& calo, KKTRK& kktrk, CCHandle cchandle, KKCALOHITCOL& hits) const { - double crystalLength = calo.caloInfo().getDouble("crystalZLength"); + double crystalLength = calo.G4Info().get("crystalZLength"); auto const& ptraj = kktrk.fitTraj(); auto cccol = cchandle.product(); double edep(-1.0); @@ -588,9 +588,9 @@ namespace mu2e { // loop over disks to decide which are worth testing std::array test{false,false}; for(unsigned idisk=0; idisk < 2; ++idisk){ - auto ffpos = calo.geomUtil().mu2eToTracker(calo.disk(idisk).geomInfo().frontFaceCenter()); - double rmin = calo.disk(idisk).geomInfo().innerEnvelopeR() - maxCaloDoca_; - double rmax = calo.disk(idisk).geomInfo().outerEnvelopeR() + maxCaloDoca_; + auto ffpos = calo.mu2eToTracker(calo.disk(idisk).diskInfo().frontFaceCenter()); + double rmin = calo.disk(idisk).diskInfo().innerEnvelopeR() - maxCaloDoca_; + double rmax = calo.disk(idisk).diskInfo().outerEnvelopeR() + maxCaloDoca_; // test at both faces; if the track is in the right area, test the clusters on this disk // Replace this with an intersection with the calo face TODO for(int iface=0; iface<2; ++iface){ @@ -740,7 +740,7 @@ namespace mu2e { // calculate the unbiased time residual Residual ctres = calohit->residual(0); // calculate the cluster depth = distance along the crystal axis from the POCA to the back face of this disk (where the SiPM sits) - double backz = calo.geomUtil().mu2eToTracker(calo.disk(calohit->caloCluster()->diskID()).geomInfo().backFaceCenter()).z(); + double backz = calo.mu2eToTracker(calo.disk(calohit->caloCluster()->diskID()).diskInfo().backFaceCenter()).z(); // calculate the distance from POCA to the SiPM, along the crystal (Z) direction, and projected along the track float clen = backz-ca.sensorPoca().Z(); float trklen = clen/ca.particleTraj().direction(ca.particleToca()).Z(); diff --git a/TrkPatRec/src/RobustMultiHelixFinder_module.cc b/TrkPatRec/src/RobustMultiHelixFinder_module.cc index 8e3ca8deb0..9ffcfb68b0 100644 --- a/TrkPatRec/src/RobustMultiHelixFinder_module.cc +++ b/TrkPatRec/src/RobustMultiHelixFinder_module.cc @@ -614,8 +614,8 @@ chi2dXY = bestHelix.fita_zt_; float caloPhi(-999),caloZ(-9999); if (circle.caloPtr_ && circle.caloPtr_->energyDep()>ccMinEnergy_){ const auto& caloCluster = *(circle.caloPtr_); - const auto caloPosInMu2e = cal_->geomUtil().diskFFToMu2e(caloCluster.diskID(),caloCluster.cog3Vector()); - const auto caloPosInTrk = cal_->geomUtil().mu2eToTracker(caloPosInMu2e); + const auto caloPosInMu2e = cal_->diskFFToMu2e(caloCluster.diskID(),caloCluster.cog3Vector()); + const auto caloPosInTrk = cal_->mu2eToTracker(caloPosInMu2e); caloZ = caloPosInTrk.z(); caloPhi = polyAtan2(caloPosInTrk.y()-circle.y_,caloPosInTrk.x()-circle.x_); } @@ -684,8 +684,8 @@ chi2dXY = bestHelix.fita_zt_; if (circle.caloPtr_ && circle.caloPtr_->energyDep()>ccMinEnergy_){ const auto& caloCluster = *(circle.caloPtr_); - const auto caloPosInMu2e = cal_->geomUtil().diskFFToMu2e(caloCluster.diskID(),caloCluster.cog3Vector()); - const auto caloPosInTrk = cal_->geomUtil().mu2eToTracker(caloPosInMu2e); + const auto caloPosInMu2e = cal_->diskFFToMu2e(caloCluster.diskID(),caloCluster.cog3Vector()); + const auto caloPosInTrk = cal_->mu2eToTracker(caloPosInMu2e); float caloPhi = polyAtan2(caloPosInTrk.y()-circle.y_,caloPosInTrk.x()-circle.x_); int nloop = round((caloPhi-(caloPosInTrk.z()-circle.fitb_zp_)/circle.fita_zp_)/6.283185); float phiLoop = caloPhi - nloop*6.283185; @@ -712,8 +712,8 @@ chi2dXY = bestHelix.fita_zt_; if (circle.caloPtr_ && circle.caloPtr_->energyDep()>ccMinEnergy_){ const auto& caloCluster = *(circle.caloPtr_); - const auto posInMu2e = cal_->geomUtil().diskFFToMu2e(caloCluster.diskID(),caloCluster.cog3Vector()); - const auto caloClusterPos = cal_->geomUtil().mu2eToTracker(posInMu2e); + const auto posInMu2e = cal_->diskFFToMu2e(caloCluster.diskID(),caloCluster.cog3Vector()); + const auto caloClusterPos = cal_->mu2eToTracker(posInMu2e); tzFitter.add(caloClusterPos.z(),caloCluster.time(),ccWeight_); } diff --git a/TrkReco/src/RobustHelixFit.cc b/TrkReco/src/RobustHelixFit.cc index 95f865ef51..0d2e1e09c0 100644 --- a/TrkReco/src/RobustHelixFit.cc +++ b/TrkReco/src/RobustHelixFit.cc @@ -1273,7 +1273,7 @@ void RobustHelixFit::findAGE(RobustHelixFinderData const& HelixData, XYZVectorF { // mu2e::GeomHandle ch; // const Calorimeter* calo = ch.get(); - XYZVectorF cog(_calorimeter->geomUtil().mu2eToTracker(_calorimeter->geomUtil().diskFFToMu2e(HelixData._hseed.caloCluster()->diskID(),HelixData._hseed.caloCluster()->cog3Vector()))); + XYZVectorF cog(_calorimeter->mu2eToTracker(_calorimeter->diskFFToMu2e(HelixData._hseed.caloCluster()->diskID(),HelixData._hseed.caloCluster()->cog3Vector()))); float rad = sqrtf(XYZVectorF(cog - center).perp2()); radii.push_back(make_pair(rad,_ccwt)); } diff --git a/TrkReco/src/TrkTimeCalculator.cc b/TrkReco/src/TrkTimeCalculator.cc index 98e1201d4b..b2a7b33940 100644 --- a/TrkReco/src/TrkTimeCalculator.cc +++ b/TrkReco/src/TrkTimeCalculator.cc @@ -50,7 +50,7 @@ namespace mu2e double TrkTimeCalculator::caloClusterTime(CaloCluster const& cc,double pitch) const { mu2e::GeomHandle ch; - CLHEP::Hep3Vector cog = ch->geomUtil().mu2eToTracker(ch->geomUtil().diskToMu2e( cc.diskID(), cc.cog3Vector())); + CLHEP::Hep3Vector cog = ch->mu2eToTracker(ch->diskToMu2e( cc.diskID(), cc.cog3Vector())); return cc.time() - timeOfFlightTimeOffset(cog.z()+_caloZOffset,pitch) + trkToCaloTimeOffset(); } From 3cba7b165f73a10303a94a32f018fd0fa4ba971d Mon Sep 17 00:00:00 2001 From: echenard Date: Wed, 29 Jul 2026 14:58:34 -0500 Subject: [PATCH 2/3] Additional CMake + overlap fix --- CalorimeterGeom/CMakeLists.txt | 1 - Mu2eG4/src/constructDiskCalorimeter.cc | 2 +- 2 files changed, 1 insertion(+), 2 deletions(-) diff --git a/CalorimeterGeom/CMakeLists.txt b/CalorimeterGeom/CMakeLists.txt index 25b8571ac2..e37eb1974b 100644 --- a/CalorimeterGeom/CMakeLists.txt +++ b/CalorimeterGeom/CMakeLists.txt @@ -1,6 +1,5 @@ cet_make_library( SOURCE - src/CaloUtil.cc src/DiskCalorimeter.cc src/Disk.cc src/DiskInfo.cc diff --git a/Mu2eG4/src/constructDiskCalorimeter.cc b/Mu2eG4/src/constructDiskCalorimeter.cc index 13a1485332..2ece9c83e7 100644 --- a/Mu2eG4/src/constructDiskCalorimeter.cc +++ b/Mu2eG4/src/constructDiskCalorimeter.cc @@ -480,7 +480,7 @@ namespace mu2e { const double crystalDZ = cal.G4Info().get("crystalZLength")/2.0; const double crystalCapDZ = cal.G4Info().get("crystalCapZLength")/2.0; const double wrapperHalfThick = cal.G4Info().get("wrapperThickness")/2.0; - const double wrapperDXY = crystalDXY + 2*wrapperHalfThick; + const double wrapperDXY = crystalDXY + 2*wrapperHalfThick - 1e-5; const double wrapperDZ = crystalDZ+crystalCapDZ; const std::vector caphriCystalId = cal.G4Info().get>("caphriCrystalId"); From 35d9c05d893f903475418c6a85ac347c2d81c8e6 Mon Sep 17 00:00:00 2001 From: echenard Date: Thu, 30 Jul 2026 12:06:37 -0500 Subject: [PATCH 3/3] Fixed bug + minor cleanup --- CaloMC/fcl/prolog.fcl | 1 - CaloMC/src/CaloShowerStepMaker_module.cc | 3 --- CalorimeterGeom/CMakeLists.txt | 1 + CalorimeterGeom/src/Disk.cc | 2 +- CalorimeterGeom/src/DiskCalorimeter.cc | 9 ++++++--- GeometryService/src/DiskCalorimeterMaker.cc | 11 +++++++---- 6 files changed, 15 insertions(+), 12 deletions(-) diff --git a/CaloMC/fcl/prolog.fcl b/CaloMC/fcl/prolog.fcl index 35b46ec002..8db3b1935d 100644 --- a/CaloMC/fcl/prolog.fcl +++ b/CaloMC/fcl/prolog.fcl @@ -38,7 +38,6 @@ CaloMC : { @table::CaloMC numZSlices : 20 deltaTime : 0.2 caloStepPointCollection : ["g4run:calorimeter"] - usePhysVolInfo : false physVolInfoInput : "g4run" caloMaterial : ["G4_CESIUM_IODIDE", "Polyethylene092","CarbonFiber"] compressData : true diff --git a/CaloMC/src/CaloShowerStepMaker_module.cc b/CaloMC/src/CaloShowerStepMaker_module.cc index be5513302e..98246c9b2f 100644 --- a/CaloMC/src/CaloShowerStepMaker_module.cc +++ b/CaloMC/src/CaloShowerStepMaker_module.cc @@ -97,7 +97,6 @@ namespace mu2e { fhicl::Atom physVolInfoInput { Name("physVolInfoInput"), Comment("Physics volume token names") }; fhicl::Atom numZSlices { Name("numZSlices"), Comment("Number of crystal longitudinal slices ") }; fhicl::Atom deltaTime { Name("deltaTime"), Comment("Max time difference to be inside a ShowerStep") }; - fhicl::Atom usePhysVolInfo { Name("usePhysVolInfo"), Comment("Use Physical Info volume names") }; fhicl::Sequence caloMaterial { Name("caloMaterial"), Comment("List of calo material names") }; fhicl::Atom compressData { Name("compressData"), Comment("Compress stepPointMC and SimParticles in crystal") }; fhicl::Atom eDepThreshold { Name("eDepThreshold"), Comment("Threshold on energy deposited by SimParticle to keep it") }; @@ -128,7 +127,6 @@ namespace mu2e { std::vector calorimeterStepPoints_; art::InputTag physVolInfoInput_; std::set mapPhysVol_; - bool usePhysVol_; std::vector caloMaterial_; int numZSlices_; double deltaTime_; @@ -156,7 +154,6 @@ namespace mu2e { art::EDProducer{config}, calorimeterStepPoints_(config().caloStepPointCollection()), physVolInfoInput_ (config().physVolInfoInput()), - usePhysVol_ (config().usePhysVolInfo()), caloMaterial_ (config().caloMaterial()), numZSlices_ (config().numZSlices()), deltaTime_ (config().deltaTime()), diff --git a/CalorimeterGeom/CMakeLists.txt b/CalorimeterGeom/CMakeLists.txt index e37eb1974b..1a3c82600a 100644 --- a/CalorimeterGeom/CMakeLists.txt +++ b/CalorimeterGeom/CMakeLists.txt @@ -14,6 +14,7 @@ cet_make_library( CLHEP::CLHEP cetlib_except::cetlib_except ) +configure_file(${CMAKE_CURRENT_SOURCE_DIR}/data/crystalPos.txt ${CURRENT_BINARY_DIR} data/crystalPos.txt) install_source(SUBDIRS src) install_headers(USE_PROJECT_NAME SUBDIRS inc) diff --git a/CalorimeterGeom/src/Disk.cc b/CalorimeterGeom/src/Disk.cc index 1cbd5e9b32..a64aad1a5f 100644 --- a/CalorimeterGeom/src/Disk.cc +++ b/CalorimeterGeom/src/Disk.cc @@ -1,5 +1,5 @@ // -// Create a disk and fills it with crystals. +// Create a disk and fills it with crystals. Crystal numbering is all local here // // Original author B Echenard // diff --git a/CalorimeterGeom/src/DiskCalorimeter.cc b/CalorimeterGeom/src/DiskCalorimeter.cc index d67aacd9be..355bedc753 100644 --- a/CalorimeterGeom/src/DiskCalorimeter.cc +++ b/CalorimeterGeom/src/DiskCalorimeter.cc @@ -14,7 +14,8 @@ namespace mu2e { ProditionsEntity(Calorimeter::cxname), disks_(), crystals_(), - G4Info_() + G4Info_(), + trackerCenter_() { rebuildCrystalPtrs(); } @@ -22,7 +23,8 @@ namespace mu2e { DiskCalorimeter::DiskCalorimeter(const DiskCalorimeter& rhs) : ProditionsEntity(rhs), disks_(rhs.disks_), - G4Info_(rhs.G4Info_) + G4Info_(rhs.G4Info_), + trackerCenter_(rhs.trackerCenter_) { rebuildCrystalPtrs(); } @@ -30,7 +32,8 @@ namespace mu2e { DiskCalorimeter::DiskCalorimeter(DiskCalorimeter&& rhs) noexcept : ProditionsEntity(std::move(rhs)), disks_(std::move(rhs.disks_)), - G4Info_(std::move(rhs.G4Info_)) + G4Info_(std::move(rhs.G4Info_)), + trackerCenter_(std::move(rhs.trackerCenter_)) { rebuildCrystalPtrs(); } diff --git a/GeometryService/src/DiskCalorimeterMaker.cc b/GeometryService/src/DiskCalorimeterMaker.cc index 6032890336..eb2f2f337d 100644 --- a/GeometryService/src/DiskCalorimeterMaker.cc +++ b/GeometryService/src/DiskCalorimeterMaker.cc @@ -58,7 +58,6 @@ namespace mu2e { calo_->G4Info_.set("caloFEBRadiusOut", config.getDouble("calorimeter.caloFEBRadiusOut") ); //outer radius of FEB enveloppe calo_->G4Info_.set("caloMotherZ0", config.getDouble("calorimeter.caloMotherZ0") ); //upstream Z of calorimeter enveloppe volume calo_->G4Info_.set("caloMotherZ1", config.getDouble("calorimeter.caloMotherZ1") ); //downstream Z of calorimeter enveloppe volume - calo_->G4Info_.set("caloMotherXorig", config.getDouble("mu2e.solenoidOffset") ); //x coordinate center of calo mother volume calo_->G4Info_.set("vdThickness", config.getDouble("calorimeter.vdThickness") ); //virtual detector thickness calo_->G4Info_.set("diskInAlRingRIn", config.getDouble("calorimeter.diskInAlRingRIn") ); //inner radius of Al ring inside disk @@ -255,10 +254,14 @@ namespace mu2e { Crystal& thisCrystal = thisDisk.crystal(icry); calo_->crystals_.push_back(&thisCrystal); - //precompute the neighbors in the global frame + //set ID and pre-compute neighbors - all need to be in global frame thisCrystal.setID(icry+crystalOffset); - thisCrystal.setNeighbors(thisDisk.neighbors(icry,1)); - thisCrystal.setNextNeighbors(thisDisk.neighbors(icry,2)); + auto neighbors = thisDisk.neighbors(icry,1); //local + auto nextNeighbors = thisDisk.neighbors(icry,2); //local + for (auto& n : neighbors){ n+=crystalOffset;} //global + for (auto& n : nextNeighbors){ n+=crystalOffset;} //global + thisCrystal.setNeighbors(neighbors); //this is local + thisCrystal.setNextNeighbors(nextNeighbors);//this is local //pre-compute the crystal position in the mu2e (global) frame CLHEP::Hep3Vector globalPosition = thisDisk.diskInfo().origin() + thisDisk.diskInfo().inverseRotation()*(thisCrystal.localPosition());