/data/doxygen/doxygen-1.7.3/gen/CMSSW_4_2_8/src/TrackingTools/TrackAssociator/src/TrackDetMatchInfo.cc

Go to the documentation of this file.

#include <map>
#include "TrackingTools/TrackAssociator/interface/TrackDetMatchInfo.h"
#include "TrackingTools/TrackAssociator/interface/DetIdInfo.h"
#include "DataFormats/EcalDetId/interface/EBDetId.h"
#include "DataFormats/EcalDetId/interface/EEDetId.h"
#include "DataFormats/CaloTowers/interface/CaloTowerDetId.h"
#include "Geometry/CaloGeometry/interface/CaloSubdetectorGeometry.h"
#include "Geometry/CaloGeometry/interface/CaloCellGeometry.h"
#include "FWCore/MessageLogger/interface/MessageLogger.h"
#include "DataFormats/Math/interface/LorentzVector.h"
#include "Math/VectorUtil.h"
#include <algorithm>



std::string TrackDetMatchInfo::dumpGeometry( const DetId& id )
{
   if ( ! caloGeometry.isValid() || 
        ! caloGeometry->getSubdetectorGeometry(id) ||
        ! caloGeometry->getSubdetectorGeometry(id)->getGeometry(id) ) {
      throw cms::Exception("FatalError")  << "Failed to access geometry for DetId: " << id.rawId();
   }
   std::ostringstream oss;

   const CaloCellGeometry::CornersVec& points = caloGeometry->getSubdetectorGeometry(id)->getGeometry(id)->getCorners();
   for( CaloCellGeometry::CornersVec::const_iterator point = points.begin();
       point != points.end(); ++point)
     oss << "(" << point->z() << ", " << point->perp() << ", " << point->eta() << ", " << point->phi() << "), \t";
   return oss.str();
}


GlobalPoint TrackDetMatchInfo::getPosition( const DetId& id)
{
   // this part might be slow
   if ( ! caloGeometry.isValid() || 
        ! caloGeometry->getSubdetectorGeometry(id) ||
        ! caloGeometry->getSubdetectorGeometry(id)->getGeometry(id) ) {
      throw cms::Exception("FatalError") << "Failed to access geometry for DetId: " << id.rawId();
      return GlobalPoint(0,0,0);
   }
   return caloGeometry->getSubdetectorGeometry(id)->getGeometry(id)->getPosition();
}


double TrackDetMatchInfo::crossedEnergy( EnergyType type )
{
   double energy(0);
   switch (type) {
    case EcalRecHits:
        {
           for(std::vector<const EcalRecHit*>::const_iterator hit=crossedEcalRecHits.begin(); hit!=crossedEcalRecHits.end(); hit++)
             energy += (*hit)->energy();
        }
      break;
    case HcalRecHits:
        {
           for(std::vector<const HBHERecHit*>::const_iterator hit = crossedHcalRecHits.begin(); hit != crossedHcalRecHits.end(); hit++)
             energy += (*hit)->energy();
        }
      break;
    case HORecHits:
        {
           for(std::vector<const HORecHit*>::const_iterator hit = crossedHORecHits.begin(); hit != crossedHORecHits.end(); hit++)
             energy += (*hit)->energy();
        }
      break;
    case TowerTotal:
        {
           for(std::vector<const CaloTower*>::const_iterator hit=crossedTowers.begin(); hit!=crossedTowers.end(); hit++)
             energy += (*hit)->energy();
        }
      break;
    case TowerEcal:
        {
           for(std::vector<const CaloTower*>::const_iterator hit=crossedTowers.begin(); hit!=crossedTowers.end(); hit++)
             energy += (*hit)->emEnergy();
        }
      break;
    case TowerHcal:
        {
           for(std::vector<const CaloTower*>::const_iterator hit=crossedTowers.begin(); hit!=crossedTowers.end(); hit++)
             energy += (*hit)->hadEnergy();
        }
      break;
    case TowerHO:
        {
           for(std::vector<const CaloTower*>::const_iterator hit=crossedTowers.begin(); hit!=crossedTowers.end(); hit++)
             energy += (*hit)->outerEnergy();
        }
      break;
    default:
      throw cms::Exception("FatalError") << "Unknown calo energy type: " << type;
   }
   return energy;
}

bool TrackDetMatchInfo::insideCone(const DetId& id, const double dR) {
   GlobalPoint idPosition = getPosition(id);
   if (idPosition.mag()<0.01) return false;
   
   math::XYZVector idPositionRoot( idPosition.x(), idPosition.y(), idPosition.z() );
   math::XYZVector trackP3( stateAtIP.momentum().x(), stateAtIP.momentum().y(), stateAtIP.momentum().z() );
   return ROOT::Math::VectorUtil::DeltaR(trackP3, idPositionRoot) < 0.5;
}

double TrackDetMatchInfo::coneEnergy( double dR, EnergyType type )
{
   double energy(0);
   switch (type) {
    case EcalRecHits:
        {
           for(std::vector<const EcalRecHit*>::const_iterator hit=ecalRecHits.begin(); hit!=ecalRecHits.end(); hit++)
             if (insideCone((*hit)->detid(),dR)) energy += (*hit)->energy();
        }
      break;
    case HcalRecHits:
        {
           for(std::vector<const HBHERecHit*>::const_iterator hit = hcalRecHits.begin(); hit != hcalRecHits.end(); hit++)
             if (insideCone((*hit)->detid(),dR)) energy += (*hit)->energy();
        }
      break;
    case HORecHits:
        {
           for(std::vector<const HORecHit*>::const_iterator hit = hoRecHits.begin(); hit != hoRecHits.end(); hit++)
             if (insideCone((*hit)->detid(),dR)) energy += (*hit)->energy();
        }
      break;
    case TowerTotal:
        {
           for(std::vector<const CaloTower*>::const_iterator hit=crossedTowers.begin(); hit!=crossedTowers.end(); hit++)
             if (insideCone((*hit)->id(),dR)) energy += (*hit)->energy();
        }
      break;
    case TowerEcal:
        {
           for(std::vector<const CaloTower*>::const_iterator hit=crossedTowers.begin(); hit!=crossedTowers.end(); hit++)
             if (insideCone((*hit)->id(),dR)) energy += (*hit)->emEnergy();
        }
      break;
    case TowerHcal:
        {
           for(std::vector<const CaloTower*>::const_iterator hit=crossedTowers.begin(); hit!=crossedTowers.end(); hit++)
             if (insideCone((*hit)->id(),dR)) energy += (*hit)->hadEnergy();
        }
      break;
    case TowerHO:
        {
           for(std::vector<const CaloTower*>::const_iterator hit=crossedTowers.begin(); hit!=crossedTowers.end(); hit++)
             if (insideCone((*hit)->id(),dR)) energy += (*hit)->outerEnergy();
        }
      break;
    default:
      throw cms::Exception("FatalError") << "Unknown calo energy type: " << type;
   }
   return energy;
}



double TrackDetMatchInfo::nXnEnergy(const DetId& id, EnergyType type, int gridSize)
{
   double energy(0);
   if ( id.rawId() == 0 ) return 0.;
   switch (type)  {
    case TowerTotal:
    case TowerHcal:
    case TowerEcal:
    case TowerHO:
        {
           if ( id.det() != DetId::Calo ) {
              throw cms::Exception("FatalError") << "Wrong DetId. Expected CaloTower, but found:\n" <<
                DetIdInfo::info(id)<<"\n";
           }
           CaloTowerDetId centerId(id);
           for(std::vector<const CaloTower*>::const_iterator hit=towers.begin(); hit!=towers.end(); hit++) {
              CaloTowerDetId neighborId((*hit)->id());
              int dEta = abs( (centerId.ieta()<0?centerId.ieta()+1:centerId.ieta() )
                              -(neighborId.ieta()<0?neighborId.ieta()+1:neighborId.ieta() ) ) ;
              int dPhi = abs( centerId.iphi()-neighborId.iphi() );
              if ( abs(72-dPhi) < dPhi ) dPhi = 72-dPhi;
              if(  dEta <= gridSize && dPhi <= gridSize ) {
                 switch (type) {
                  case TowerTotal:
                    energy += (*hit)->energy();
                    break;
                  case TowerEcal:
                    energy += (*hit)->emEnergy();
                    break;
                  case TowerHcal:
                    energy += (*hit)->hadEnergy();
                    break;
                  case TowerHO:
                    energy += (*hit)->outerEnergy();
                    break;
                  default:
                    throw cms::Exception("FatalError") << "Unknown calo tower energy type: " << type;
                 }
              }
           }
        }
      break;
    case EcalRecHits:
        {
           if( id.det() != DetId::Ecal || (id.subdetId() != EcalBarrel && id.subdetId() != EcalEndcap) ) {
              throw cms::Exception("FatalError") << "Wrong DetId. Expected EcalBarrel or EcalEndcap, but found:\n" <<
                DetIdInfo::info(id)<<"\n";
           }
           // Since the ECAL granularity is small and the gap between EE and EB is significant,
           // energy is computed only within the system that contains the central element
           if( id.subdetId() == EcalBarrel ) {
              EBDetId centerId(id);
              for(std::vector<const EcalRecHit*>::const_iterator hit=ecalRecHits.begin(); hit!=ecalRecHits.end(); hit++) {
                 if ((*hit)->id().subdetId() != EcalBarrel) continue;
                 EBDetId neighborId((*hit)->id());
                 int dEta = abs( (centerId.ieta()<0?centerId.ieta()+1:centerId.ieta() )
                                 -(neighborId.ieta()<0?neighborId.ieta()+1:neighborId.ieta() ) ) ;
                 int dPhi = abs( centerId.iphi()-neighborId.iphi() );
                 if ( abs(360-dPhi) < dPhi ) dPhi = 360-dPhi;
                 if(  dEta <= gridSize && dPhi <= gridSize ) {
                    energy += (*hit)->energy();
                 }
              }
           } else {
              // Endcap
              EEDetId centerId(id);
              for(std::vector<const EcalRecHit*>::const_iterator hit=ecalRecHits.begin(); hit!=ecalRecHits.end(); hit++) {
                 if ((*hit)->id().subdetId() != EcalEndcap) continue;
                 EEDetId neighborId((*hit)->id());
                 if(  centerId.zside() == neighborId.zside() && 
                      abs(centerId.ix()-neighborId.ix()) <= gridSize && 
                      abs(centerId.iy()-neighborId.iy()) <= gridSize ) {
                    energy += (*hit)->energy();
                 }
              }
           }
        }
      break;
    case HcalRecHits:
        {
           if( id.det() != DetId::Hcal || (id.subdetId() != HcalBarrel && id.subdetId() != HcalEndcap) ) {
              throw cms::Exception("FatalError") << "Wrong DetId. Expected HE or HB, but found:\n" <<
                DetIdInfo::info(id)<<"\n";
           }
           HcalDetId centerId(id);
           for(std::vector<const HBHERecHit*>::const_iterator hit=hcalRecHits.begin(); hit!=hcalRecHits.end(); hit++) {
              HcalDetId neighborId((*hit)->id());
              int dEta = abs( (centerId.ieta()<0?centerId.ieta()+1:centerId.ieta() )
                              -(neighborId.ieta()<0?neighborId.ieta()+1:neighborId.ieta() ) ) ;
              int dPhi = abs( centerId.iphi()-neighborId.iphi() );
              if ( abs(72-dPhi) < dPhi ) dPhi = 72-dPhi;
              if(  dEta <= gridSize && dPhi <= gridSize ){
                 energy += (*hit)->energy();
              }
           }
        }
      break;
    case HORecHits:
        {
           if( id.det() != DetId::Hcal || (id.subdetId() != HcalOuter) ) {
              throw cms::Exception("FatalError") << "Wrong DetId. Expected HO, but found:\n" <<
                DetIdInfo::info(id)<<"\n";
           }
           HcalDetId centerId(id);
           for(std::vector<const HORecHit*>::const_iterator hit=hoRecHits.begin(); hit!=hoRecHits.end(); hit++) {
              HcalDetId neighborId((*hit)->id());
              int dEta = abs( (centerId.ieta()<0?centerId.ieta()+1:centerId.ieta() )
                              -(neighborId.ieta()<0?neighborId.ieta()+1:neighborId.ieta() ) ) ;
              int dPhi = abs( centerId.iphi()-neighborId.iphi() );
              if ( abs(72-dPhi) < dPhi ) dPhi = 72-dPhi;
              if(  dEta <= gridSize && dPhi <= gridSize ) {
                 energy += (*hit)->energy();
              }
           }
        }
      break;
    default:
      throw cms::Exception("FatalError") << "Unkown or not implemented energy type requested, type:" << type;
   }
   return energy;
}

double TrackDetMatchInfo::nXnEnergy(EnergyType type, int gridSize)
{
   switch (type)  {
    case TowerTotal:
    case TowerHcal:
    case TowerEcal:
    case TowerHO:
      if( crossedTowerIds.empty() ) return 0;
      return nXnEnergy(crossedTowerIds.front(), type, gridSize);
      break;
    case EcalRecHits:
      if( crossedEcalIds.empty() ) return 0;
      return nXnEnergy(crossedEcalIds.front(), type, gridSize);
      break;
    case HcalRecHits:
      if( crossedHcalIds.empty() ) return 0;
      return nXnEnergy(crossedHcalIds.front(), type, gridSize);
      break;
    case HORecHits:
      if( crossedHOIds.empty() ) return 0;
      return nXnEnergy(crossedHOIds.front(), type, gridSize);
      break;
    default:
      throw cms::Exception("FatalError") << "Unkown or not implemented energy type requested, type:" << type;
   }
   return -999;
}


TrackDetMatchInfo::TrackDetMatchInfo():
       trkGlobPosAtEcal(0,0,0)
     , trkGlobPosAtHcal(0,0,0)
     , trkGlobPosAtHO(0,0,0)
     , trkMomAtEcal(0,0,0)
     , trkMomAtHcal(0,0,0)
     , trkMomAtHO(0,0,0)
     , isGoodEcal(false)
     , isGoodHcal(false)
     , isGoodCalo(false)
     , isGoodHO(false)
     , isGoodMuon(false)
     , simTrack(0)
     , ecalTrueEnergy(-999)
     , hcalTrueEnergy(-999)
{
}

DetId TrackDetMatchInfo::findMaxDeposition( EnergyType type )
{
   DetId id;
   float maxEnergy = -9999;
   switch (type) {
    case EcalRecHits:
        {
           for(std::vector<const EcalRecHit*>::const_iterator hit=ecalRecHits.begin(); hit!=ecalRecHits.end(); hit++)
             if ( (*hit)->energy() > maxEnergy ) {
                maxEnergy = (*hit)->energy();
                id = (*hit)->detid();
             }
        }
      break;
    case HcalRecHits:
        {
           for(std::vector<const HBHERecHit*>::const_iterator hit=hcalRecHits.begin(); hit!=hcalRecHits.end(); hit++)
             if ( (*hit)->energy() > maxEnergy ) {
                maxEnergy = (*hit)->energy();
                id = (*hit)->detid();
             }
        }
      break;
    case HORecHits:
        {
           for(std::vector<const HORecHit*>::const_iterator hit=hoRecHits.begin(); hit!=hoRecHits.end(); hit++)
             if ( (*hit)->energy() > maxEnergy ) {
                maxEnergy = (*hit)->energy();
                id = (*hit)->detid();
             }
        }
      break;
    case TowerTotal:
    case TowerEcal:
    case TowerHcal:
        {
           for(std::vector<const CaloTower*>::const_iterator hit=towers.begin(); hit!=towers.end(); hit++)
             {
                double energy = 0;
                switch (type) {
                 case TowerTotal:
                   energy = (*hit)->energy();
                   break;
                 case TowerEcal:
                   energy = (*hit)->emEnergy();
                   break;
                 case TowerHcal:
                   energy = (*hit)->hadEnergy();
                   break;
                 case TowerHO:
                   energy = (*hit)->energy();
                   break;
                 default:
                   throw cms::Exception("FatalError") << "Unknown calo tower energy type: " << type;
                }
                if ( energy > maxEnergy ) {
                   maxEnergy = energy;
                   id = (*hit)->id();
                }
             }
        }
    default:
      throw cms::Exception("FatalError") << "Maximal energy deposition: unkown or not implemented energy type requested, type:" << type;
   }
   return id;
}

DetId TrackDetMatchInfo::findMaxDeposition(const DetId& id, EnergyType type, int gridSize)
{
   double energy_max(0);
   DetId id_max;
   if ( id.rawId() == 0 ) return id_max;
   switch (type)  {
    case TowerTotal:
    case TowerHcal:
    case TowerEcal:
    case TowerHO:
        {
           if ( id.det() != DetId::Calo ) {
              throw cms::Exception("FatalError") << "Wrong DetId. Expected CaloTower, but found:\n" <<
                DetIdInfo::info(id)<<"\n";
           }
           CaloTowerDetId centerId(id);
           for(std::vector<const CaloTower*>::const_iterator hit=towers.begin(); hit!=towers.end(); hit++) {
              CaloTowerDetId neighborId((*hit)->id());
              int dEta = abs( (centerId.ieta()<0?centerId.ieta()+1:centerId.ieta() )
                              -(neighborId.ieta()<0?neighborId.ieta()+1:neighborId.ieta() ) ) ;
              int dPhi = abs( centerId.iphi()-neighborId.iphi() );
              if ( abs(72-dPhi) < dPhi ) dPhi = 72-dPhi;
              if(  dEta <= gridSize && dPhi <= gridSize ) {
                 switch (type) {
                  case TowerTotal:
                    if ( energy_max < (*hit)->energy() ){
                       energy_max = (*hit)->energy();
                       id_max = (*hit)->id();
                    }
                    break;
                  case TowerEcal:
                    if ( energy_max < (*hit)->emEnergy() ){
                       energy_max = (*hit)->emEnergy();
                       id_max = (*hit)->id();
                    }
                    break;
                  case TowerHcal:
                    if ( energy_max < (*hit)->hadEnergy() ){
                       energy_max = (*hit)->hadEnergy();
                       id_max = (*hit)->id();
                    }
                    break;
                  case TowerHO:
                    if ( energy_max < (*hit)->outerEnergy() ){
                       energy_max = (*hit)->outerEnergy();
                       id_max = (*hit)->id();
                    }
                    break;
                  default:
                    throw cms::Exception("FatalError") << "Unknown calo tower energy type: " << type;
                 }
              }
           }
        }
      break;
    case EcalRecHits:
        {
           if( id.det() != DetId::Ecal || (id.subdetId() != EcalBarrel && id.subdetId() != EcalEndcap) ) {
              throw cms::Exception("FatalError") << "Wrong DetId. Expected EcalBarrel or EcalEndcap, but found:\n" <<
                DetIdInfo::info(id)<<"\n";
           }
           // Since the ECAL granularity is small and the gap between EE and EB is significant,
           // energy is computed only within the system that contains the central element
           if( id.subdetId() == EcalBarrel ) {
              EBDetId centerId(id);
              for(std::vector<const EcalRecHit*>::const_iterator hit=ecalRecHits.begin(); hit!=ecalRecHits.end(); hit++) {
                 if ((*hit)->id().subdetId() != EcalBarrel) continue;
                 EBDetId neighborId((*hit)->id());
                 int dEta = abs( (centerId.ieta()<0?centerId.ieta()+1:centerId.ieta() )
                                 -(neighborId.ieta()<0?neighborId.ieta()+1:neighborId.ieta() ) ) ;
                 int dPhi = abs( centerId.iphi()-neighborId.iphi() );
                 if ( abs(360-dPhi) < dPhi ) dPhi = 360-dPhi;
                 if(  dEta <= gridSize && dPhi <= gridSize ) {
                    if ( energy_max < (*hit)->energy() ){
                       energy_max = (*hit)->energy();
                       id_max = (*hit)->id();
                    }
                 }
              }
           } else {
              // Endcap
              EEDetId centerId(id);
              for(std::vector<const EcalRecHit*>::const_iterator hit=ecalRecHits.begin(); hit!=ecalRecHits.end(); hit++) {
                 if ((*hit)->id().subdetId() != EcalEndcap) continue;
                 EEDetId neighborId((*hit)->id());
                 if(  centerId.zside() == neighborId.zside() && 
                      abs(centerId.ix()-neighborId.ix()) <= gridSize && 
                      abs(centerId.iy()-neighborId.iy()) <= gridSize ) {
                    if ( energy_max < (*hit)->energy() ){
                       energy_max = (*hit)->energy();
                       id_max = (*hit)->id();
                    }
                 }
              }
           }
        }
      break;
    case HcalRecHits:
        {
           if( id.det() != DetId::Hcal || (id.subdetId() != HcalBarrel && id.subdetId() != HcalEndcap) ) {
              throw cms::Exception("FatalError") << "Wrong DetId. Expected HE or HB, but found:\n" <<
                DetIdInfo::info(id)<<"\n";
           }
           HcalDetId centerId(id);
           for(std::vector<const HBHERecHit*>::const_iterator hit=hcalRecHits.begin(); hit!=hcalRecHits.end(); hit++) {
              HcalDetId neighborId((*hit)->id());
              int dEta = abs( (centerId.ieta()<0?centerId.ieta()+1:centerId.ieta() )
                              -(neighborId.ieta()<0?neighborId.ieta()+1:neighborId.ieta() ) ) ;
              int dPhi = abs( centerId.iphi()-neighborId.iphi() );
              if ( abs(72-dPhi) < dPhi ) dPhi = 72-dPhi;
              if(  dEta <= gridSize && dPhi <= gridSize ){
                 if ( energy_max < (*hit)->energy() ){
                    energy_max = (*hit)->energy();
                    id_max = (*hit)->id();
                 }
              }
           }
        }
      break;
    case HORecHits:
        {
           if( id.det() != DetId::Hcal || (id.subdetId() != HcalOuter) ) {
              throw cms::Exception("FatalError") << "Wrong DetId. Expected HO, but found:\n" <<
                DetIdInfo::info(id)<<"\n";
           }
           HcalDetId centerId(id);
           for(std::vector<const HORecHit*>::const_iterator hit=hoRecHits.begin(); hit!=hoRecHits.end(); hit++) {
              HcalDetId neighborId((*hit)->id());
              int dEta = abs( (centerId.ieta()<0?centerId.ieta()+1:centerId.ieta() )
                              -(neighborId.ieta()<0?neighborId.ieta()+1:neighborId.ieta() ) ) ;
              int dPhi = abs( centerId.iphi()-neighborId.iphi() );
              if ( abs(72-dPhi) < dPhi ) dPhi = 72-dPhi;
              if(  dEta <= gridSize && dPhi <= gridSize ) {
                 if ( energy_max < (*hit)->energy() ){
                    energy_max = (*hit)->energy();
                    id_max = (*hit)->id();
                 }
              }
           }
        }
      break;
    default:
      throw cms::Exception("FatalError") << "Unkown or not implemented energy type requested, type:" << type;
   }
   return id_max;
}

DetId TrackDetMatchInfo::findMaxDeposition(EnergyType type, int gridSize)
{
   DetId id_max;
   switch (type)  {
    case TowerTotal:
    case TowerHcal:
    case TowerEcal:
    case TowerHO:
      if( crossedTowerIds.empty() ) return id_max;
      return findMaxDeposition(crossedTowerIds.front(), type, gridSize);
      break;
    case EcalRecHits:
      if( crossedEcalIds.empty() ) return id_max;
      return findMaxDeposition(crossedEcalIds.front(), type, gridSize);
      break;
    case HcalRecHits:
      if( crossedHcalIds.empty() ) return id_max;
      return findMaxDeposition(crossedHcalIds.front(), type, gridSize);
      break;
    case HORecHits:
      if( crossedHOIds.empty() ) return id_max;
      return findMaxDeposition(crossedHOIds.front(), type, gridSize);
      break;
    default:
      throw cms::Exception("FatalError") << "Unkown or not implemented energy type requested, type:" << type;
   }
   return id_max;
}

// Obsolete
//



double TrackDetMatchInfo::ecalConeEnergy()
{
   return coneEnergy (999, EcalRecHits );
}

double TrackDetMatchInfo::hcalConeEnergy()
{
   return coneEnergy (999, HcalRecHits );
}

double TrackDetMatchInfo::hoConeEnergy()
{
   return coneEnergy (999, HcalRecHits );
}

double TrackDetMatchInfo::ecalCrossedEnergy()
{
   return crossedEnergy( EcalRecHits );
}

double TrackDetMatchInfo::hcalCrossedEnergy()
{
   return crossedEnergy( HcalRecHits );
}

double TrackDetMatchInfo::hoCrossedEnergy()
{
   return crossedEnergy( HORecHits );
}
      
      
int TrackDetMatchInfo::numberOfSegments() const {
   int numSegments = 0;
   for(std::vector<TAMuonChamberMatch>::const_iterator chamber=chambers.begin(); chamber!=chambers.end(); chamber++)
     numSegments += chamber->segments.size();
   return numSegments;
}

int TrackDetMatchInfo::numberOfSegmentsInStation(int station) const {
   int numSegments = 0;
   for(std::vector<TAMuonChamberMatch>::const_iterator chamber=chambers.begin(); chamber!=chambers.end(); chamber++)
     if(chamber->station()==station) numSegments += chamber->segments.size();
   return numSegments;
}

int TrackDetMatchInfo::numberOfSegmentsInStation(int station, int detector) const {
   int numSegments = 0;
   for(std::vector<TAMuonChamberMatch>::const_iterator chamber=chambers.begin(); chamber!=chambers.end(); chamber++)
     if(chamber->station()==station&&chamber->detector()==detector) numSegments += chamber->segments.size();
   return numSegments;
}

int TrackDetMatchInfo::numberOfSegmentsInDetector(int detector) const {
   int numSegments = 0;
   for(std::vector<TAMuonChamberMatch>::const_iterator chamber=chambers.begin(); chamber!=chambers.end(); chamber++)
     if(chamber->detector()==detector) numSegments += chamber->segments.size();
   return numSegments;
}