#include <RPCSensitiveDetectorCosmics.h>

Inheritance diagram for RPCSensitiveDetectorCosmics:

Collaboration diagram for RPCSensitiveDetectorCosmics:

Public Member Functions
	RPCSensitiveDetectorCosmics (const std::string &name, const std::string &hitCollectionName, unsigned int nGasGaps)
	construction/destruction
void	Initialize (G4HCofThisEvent *) override final
	member functions
G4bool	ProcessHits (G4Step , G4TouchableHistory ) override final

Private Member Functions
	FRIEND_TEST (RPCSensitiveDetectorCosmicstest, Initialize)
	FRIEND_TEST (RPCSensitiveDetectorCosmicstest, ProcessHits)

Private Attributes
std::string	m_hitCollectionName
	member data
RPCSimHitCollection *	m_myRPCHitColl {nullptr}
AtlasG4EventUserInfo *	m_g4UserEventInfo {nullptr}
const RpcHitIdHelper *	m_muonHelper
double	m_globalTime
bool	m_isGeoModel
Amg::Vector3D	m_mom
double	m_momMag
Amg::Vector3D	m_vertex
Amg::Vector3D	m_currVertex
Amg::Vector3D	m_globH

Detailed Description

Author: Andre.nosp@m.a.Di.nosp@m.Simon.nosp@m.e@ce.nosp@m.rn.ch

Class methods and properties

The method RPCSensitiveDetectorCosmics::ProcessHits is executed by the G4 kernel each time a particle crosses one of the RPC gas gaps. Navigating with the touchableHistory method GetHistoryDepth() through the hierarchy of volumes crossed by the particle, the Sensitive Detector determines the correct set of Simulation Identifiers to associate to each hit. The RPC SimIDs are 32-bit unsigned integers, built using the MuonSimEvent/RpcHitIdHelper class which inherits from the MuonHitIdHelper base class.

We describe here how each field of the identifier is determined.

1) stationName, stationEta, stationPhi: when a volume is found in the hierarchy whose name contains the substring "station", the stationName is extracted from the volume's name. stationEta and stationPhi values are calculated starting from the volume's copy number.

2) doubletR: when a volume is found in the hierarchy whose name contains the substring "rpccomponent", its copy number is used to calulate the doubletR of the hit.

3) gasGap: the substring "layer" is searched through the names of the volumes in the hierarchy. When a the volume is found, its copy number is used, together with the sign of the copy number of the "rpccomponent" to decide what is the correct gasGap value.

4) doubletPhi: when the volume with the substring "gas volume" in its name is found, its copy number is enough to determine the doubletPhi of the hit if the hit is registered in a standard chamber. For special chambers (i.e. chambers with only one gas gap per layer readout by two strip panels per direction, or chambers in the ribs) some a special attribution is done using also the stationName.

5) doubletZ: "gazGap" is the required substring in the volume's name. doubletZ attribution is fairly complicated. For standard chambers it just uses the sign of the copy number of the "rpccomponent", but severa arrangements (based on the stationName and the technology name) are needed for special chambers.

6) the doubletPhi calculated as described above needs one further correction for some chambers to take into account chamber rotation before placement in the spectrometer. This is done just before creating the hit.

notes:

1) presently, chamber efficiency is assumed to be 100% at the level of the Sensitive Detector, i.e. two hits (one in eta and one in phi) are created each time the sensitive detector is created.

2) the hits created as described above contain information about their local position in the gas gap and the simulation identifier of the strip panel

3) the strip number is not assigned by the sensitive detector, and must be calculated by the digitization algorithm.

4) points 2 and 3 are due to the necessity of not introducing any dependency on the geometry description in the Sensitive Detector.

5) the present version of the Sensitive Detector produces hits which are different depending on whether hand coded G4 geometry or GeoModel is used for the geometrical description of the setup. When hand coded geometry is used, both GeoModel and the old MuonDetDescr can be used for digitization (using a proper tag of RPC_Digitization), while if GeoModel is used in the simulation, it must be used also in the digitization.

6) for each hit, the time of flight (the G4 globalTime), is recorded and associated to the hit.

7) the RPCHit object contains: the SimID, the globalTime, the hit local position and the track barcode.

Definition at line 98 of file RPCSensitiveDetectorCosmics.h.

Constructor & Destructor Documentation

◆ RPCSensitiveDetectorCosmics()

RPCSensitiveDetectorCosmics::RPCSensitiveDetectorCosmics	(	const std::string &	name,
		const std::string &	hitCollectionName,
		unsigned int	nGasGaps )

construction/destruction

Definition at line 22 of file RPCSensitiveDetectorCosmics.cxx.

  : G4VSensitiveDetector( name )
  , m_hitCollectionName( hitCollectionName )
  , m_globalTime(0.)
  , m_isGeoModel(true)
  , m_momMag(0.)
{
  m_muonHelper = RpcHitIdHelper::GetHelper(nGasGaps);
}

Member Function Documentation

◆ FRIEND_TEST() [1/2]

RPCSensitiveDetectorCosmics::FRIEND_TEST	(	RPCSensitiveDetectorCosmicstest	,
		Initialize	)

private

◆ FRIEND_TEST() [2/2]

RPCSensitiveDetectorCosmics::FRIEND_TEST	(	RPCSensitiveDetectorCosmicstest	,
		ProcessHits	)

private

◆ Initialize()

void RPCSensitiveDetectorCosmics::Initialize ( G4HCofThisEvent * )

finaloverride

member functions

Definition at line 32 of file RPCSensitiveDetectorCosmics.cxx.

{
  m_myRPCHitColl = nullptr;
  if (auto* eventInfo = AtlasG4EventUserInfo::GetEventUserInfo()) {
    m_myRPCHitColl = eventInfo->GetHitCollectionMap()->Find<RPCSimHitCollection>(m_hitCollectionName);
    m_g4UserEventInfo = eventInfo;
  }
  if (verboseLevel>1) G4cout << "Initializing SD" << G4endl;
  // FIXME this next bit probebly only needs to be done once pre job
  // rather than once per G4Event?
/*  DetectorGeometryHelper DGHelp;
  if(  DGHelp.GeometryType("Muon") == GeoModel ){
    m_isGeoModel = true;
    if (verboseLevel>1) G4cout << "Muon Geometry is from GeoModel" << G4endl;
  } else {
    m_isGeoModel = false;
    if (verboseLevel>1) G4cout << "Muon Geometry is from pure G4" << G4endl;
  }*/
  m_mom = Amg::Vector3D(0.,0.,0.);
  m_globH = Amg::Vector3D(0.,0.,0.);
}

◆ ProcessHits()

G4bool RPCSensitiveDetectorCosmics::ProcessHits	(	G4Step *	aStep,
		G4TouchableHistory *	)

finaloverride

Definition at line 54 of file RPCSensitiveDetectorCosmics.cxx.

                                                                                 {
 
  if (!m_myRPCHitColl) {
    G4Exception("RPCSensitiveDetectorCosmics::ProcessHits", "RPCCosmicHitCollectionMissing", FatalException,
                "Hit collection not initialized; did SetupEvent run?");
    return false;
  }
  G4Track* currentTrack = aStep->GetTrack();
 
  if (currentTrack->GetDefinition()->GetPDGCharge() == 0.0) {
    if (currentTrack->GetDefinition()!=G4Geantino::GeantinoDefinition()) return true;
    // else if (currentTrack->GetDefinition()==G4ChargedGeantino::ChargedGeantinoDefinition()) return true;
  }
 
  const G4TouchableHistory* touchHist = static_cast<const G4TouchableHistory*>(aStep->GetPreStepPoint()->GetTouchable());
  G4ThreeVector           position  = aStep->GetPreStepPoint()->GetPosition();
  G4ThreeVector       postPosition  = aStep->GetPostStepPoint()->GetPosition();
  const G4AffineTransform trans     = currentTrack->GetTouchable()->GetHistory()->GetTopTransform(); // from global to local
 
  // necessary to assign correct identifiers
  int rpcIsRotated      = 0;
  // int layer             = 0;
 
  std::string tech;
  // fields for the RPC identifier construction
  std::string stationName;
  int         stationEta= 0;
  int         stationPhi= 0;
  int         doubletR= 0;
  int         doubletZ= 0;
  int         doubletPhi= 0;
  int         gasGap= 0;
  // int         eta_inversed=0;
  // int         phi_inversed=0;
  int         station_rotated=0; // tells us if the station was rotated before being positioned.
                                 // if =1 we have to correct some IDs
  int         zNeg_original=0; // tells if the station at z<0 was obtained duplicating a station at z>0
 
  // RPC hit information
  double     globalTime    = aStep->GetPreStepPoint()->GetGlobalTime();
  Amg::Vector3D localPosition = Amg::Hep3VectorToEigen( trans.TransformPoint(position) );
  Amg::Vector3D localPostPosition = Amg::Hep3VectorToEigen( trans.TransformPoint(postPosition) );
 
  // global coordinates
  G4ThreeVector globVrtx = aStep->GetPreStepPoint()->GetPosition();
 
  // distance of the hit from (0,0,0) vertex - calculated from the PreStepPoint (approximation)
  // double dist = globVrtx.mag();
  double lightspeed = 299.792458; /* in vacuo speed of light [mm/ns] */
  // double tOrigin = dist / lightspeed;
 
  G4int trackid = aStep->GetTrack()->GetTrackID();
  m_currVertex = Amg::Hep3VectorToEigen( aStep->GetTrack()->GetVertexPosition() );
 
  // for cosmics: only primary muon tracks - track momentum when first entering the spectrometer (one muon per event)
  if ((m_currVertex != m_vertex) && (trackid == 1)) {
    // after calculationg the momentum magnidude, normalize it
    m_mom = Amg::Hep3VectorToEigen( currentTrack->GetMomentum() );
    m_momMag = m_mom.mag();
    m_mom.normalize();
    // the direction of the primary mu is used to calculate the t0, the position ot the t0, m_globH, is ONE for a track
    Amg::Vector3D globVrtxFix = Amg::Hep3VectorToEigen( globVrtx );
    double AlphaGlobal = -1*(globVrtxFix[0]*m_mom[0] + globVrtxFix[1]*m_mom[1] + globVrtxFix[2]*m_mom[2])/(m_mom[0]*m_mom[0] + m_mom[1]*m_mom[1] + m_mom[2]*m_mom[2]);
    m_globH = globVrtxFix + AlphaGlobal*m_mom;
  }
  double globalDist = sqrt((m_globH[0] - globVrtx[0])*(m_globH[0] - globVrtx[0]) +
                           (m_globH[1] - globVrtx[1])*(m_globH[1] - globVrtx[1]) +
                           (m_globH[2] - globVrtx[2])*(m_globH[2] - globVrtx[2]));
  double tof = globalDist / lightspeed;
 
  int mydbZ=0;
  int mydbPMod=0;
  int mydbP=0;
  bool isAssembly = false;
  // scan geometry tree to identify the hit channel
  for (int i=touchHist->GetHistoryDepth();i>=0;i--) {
 
    std::string::size_type npos;
    std::string volName = touchHist->GetVolume(i)->GetName();
    std::string num=volName.substr(3,2);
    if(num[0]==' ') num[0]=0;
 
    // check if this station is an assembly
    if ((npos = volName.find("av_")) != std::string::npos &&
        (npos = volName.find("impr_")) != std::string::npos)  isAssembly = true;
 
    // stationName, stationEta, stationPhi
    if ((npos = volName.find("station")) != std::string::npos  && (!isAssembly)) {
 
      stationName = volName.substr(0,npos-1);
 
      int volCopyNo = touchHist->GetVolume(i)->GetCopyNo();
 
      if(abs(volCopyNo/1000)==1){
        zNeg_original=1;
        volCopyNo=volCopyNo%1000;
      }
 
      stationEta = volCopyNo/100;
      stationPhi = abs(volCopyNo%100);
 
      if(stationEta<0&&!zNeg_original) station_rotated=1;
 
      // stationName, stationEta, stationPhi
    } else if ((npos = volName.find("RPC")) != std::string::npos && isAssembly) {
      // vol name for Assembly components are
      // av_WWW_impr_XXX_Muon::BMSxMDTxx_pv_ZZZ_NAME
      //         where WWW is ass. istance nr.
      //               XXX is comp. imprint nr.
      //               BMSxMDTxx is the name of the comp. log.Vol.
      //                  x station sub-type; xx technology subtype
      //               ZZZ is the comp. number of order
      //               NAME is the comp. tag (geoIdentifierTag)
      //                    for RPCs it is SwapdbPdbRdbZ[aaa]
      //                    with aaa = doubletZ + doubletR*100 + dbphi*1000;
      //                         aaa = -aaa if iswap == -1
      //                         dbphi = 1 always but for S1 or S2 at doubletPhi = 2
      // copy numbers for Ass.components are =
      //                CopyNoBase(= geoIdentifierTag of the assembly) + geoIdentifierTag of the component
      //                geoIdentifierTag of the component = Job
      //                geoIdentifierTag of the assembly = (sideC*10000 +
      //                             mirsign*1000 + abs(zi)*100 + fi+1)*100000;
      //                             mirsign = 1 if zi<0 and !mirrored; otherwise 0
      std::string::size_type loc1,loc2;
      if ((loc1 = volName.find("Muon::")) != std::string::npos) {
        stationName = volName.substr(loc1+6,4); //type and subtype
      }
 
      int volCopyNo = touchHist->GetVolume(i)->GetCopyNo();
      int copyNrBase = int(volCopyNo/100000);
      int sideC  = int(copyNrBase/10000);
      int zi     = int((copyNrBase%1000)/100);
      int mirfl  = int((copyNrBase%10000)/1000);
      int fi     = int(copyNrBase%100);
      if (sideC == 1) zi = -zi;
      stationEta    = zi;
      stationPhi    = fi;
      zNeg_original = mirfl;
 
      if(stationEta<0&&!zNeg_original) station_rotated=1;
      tech=volName.substr(npos,5);
 
      // now get the geoIdentifierTag of the rpc components
      int gmID = 0;
      if ((loc1 = volName.find('[')) != std::string::npos) {
        if ((loc2 = volName.find(']', loc1+1)) != std::string::npos) {
          std::istringstream istrvar(volName.substr(loc1+1,loc2-loc1-1));
          istrvar>>gmID;
        }
      }
      int kk=gmID;
 
      if (kk < 0) rpcIsRotated=1;
 
      doubletR =(abs(kk)%1000)/100;
      mydbZ    = abs(int(kk%10));
      mydbPMod = abs(int(kk/1000));
      // doubleR, doubletZ
    } else if ((npos = volName.find("rpccomponent")) != std::string::npos && (!isAssembly)) {
 
      tech=volName.substr(npos-5,5);
 
      std::string::size_type loc1,loc2;
      int gmID = 0;
      if ((loc1 = volName.find('[')) != std::string::npos) {
        if ((loc2 = volName.find(']', loc1+1)) != std::string::npos) {
          std::istringstream istrvar(volName.substr(loc1+1,loc2-loc1-1));
          istrvar>>gmID;
        }
      }
      mydbZ    = abs(int(gmID%10));
      mydbPMod = abs(int(gmID/1000));
 
      int kk=touchHist->GetVolume(i)->GetCopyNo();
 
      if (kk < 0) rpcIsRotated=1;
 
      doubletR=(abs(kk)%1000)/100;
 
      // doubleR
    } else if ((npos = volName.find("rpccomponent")) != std::string::npos && (!isAssembly)) {
 
      std::string::size_type loc1,loc2;
      int gmID = 0;
      if ((loc1 = volName.find('[')) != std::string::npos) {
        if ((loc2 = volName.find(']', loc1+1)) != std::string::npos) {
          std::istringstream istrvar(volName.substr(loc1+1,loc2-loc1-1));
          istrvar>>gmID;
        }
      }
      mydbZ    = abs(int(gmID%10));
      mydbPMod = abs(int(gmID/1000));
 
      int kk=touchHist->GetVolume(i)->GetCopyNo();
 
      if (kk < 0) rpcIsRotated=1;
 
      doubletR=(abs(kk)%1000)/100;
 
    // gasGap
    } else if ((npos = volName.find("layer")) != std::string::npos) {
 
      int copyNo = touchHist->GetVolume(i)->GetCopyNo();
 
      if (copyNo == 1) {
        rpcIsRotated ? gasGap = 2 : gasGap = 1;
      } else if (copyNo ==2) {
        rpcIsRotated ? gasGap = 1 : gasGap = 2;
      }
 
      // doubletPhi
    } else if((npos = volName.find("gas volume")) != std::string::npos) {
 
      int copyNo = touchHist->GetVolume(i)->GetCopyNo();
      if (copyNo == 0) {
        doubletPhi = 1;
      } else if (copyNo == 10) {
        doubletPhi = 2;
      }
      mydbP = doubletPhi;
      int ngap_in_s=0;
      int nstrippanel_in_s=0;
      std::string::size_type loc1;
      if ((loc1 = volName.find("gg_in_s")) != std::string::npos) {
        std::istringstream istrvar(volName.substr(loc1-1,1));
        istrvar>>ngap_in_s;
      }
      if ((loc1 = volName.find("sp_in_s")) != std::string::npos) {
        std::istringstream istrvar(volName.substr(loc1-1,1));
        istrvar>>nstrippanel_in_s;
      }
      if (ngap_in_s == 1 && nstrippanel_in_s == 2){
        if(localPosition.y()>0) mydbP=2;
        else mydbP=1;
      } else if (ngap_in_s == 1 && nstrippanel_in_s == 1){
        mydbP = mydbPMod;
      }
    }
  }
  // correct wrong IDs due to station rotation. only doubletZ and doubletPhi are affected at this point
 
  if(station_rotated){
    // doubletPhi correction
    // note that this is correct also for ribs chambers
    mydbP++;
    if (mydbP>2) mydbP=1;
    // doubletZ correction: specialStations have already been correted
    // the others do not need any correction, because doubletZ increases with |Z|, thus not
    // changing in the rotation
 
    // strip numbering: if the station is rotated both eta and phi directions get inversed.
    // commented out for geomodel!
 
    if(!m_isGeoModel){
      localPosition.z() = -localPosition.z();
      localPosition.y() = -localPosition.y();
      localPostPosition.z() = -localPostPosition.z();
      localPostPosition.y() = -localPostPosition.y();
    }
  }
 
  // further correction on the eta direction due to rpc component rotation
  // commented for geomodel!
  if(!rpcIsRotated&&!m_isGeoModel) {
    localPosition.z() = -localPosition.z();
    localPostPosition.z() = -localPostPosition.z();
  }
 
  // now we have the position in the gas gap, with correct axis orientation, and the offlineIDs *of the strip panel*
  // nstrip is supposed to be calculated by RPC_Digitizer.
 
  // construct two (one in eta and one in phi) new RPC hits and store them in hit collection
 
  //construct the hit identifiers
 
  if (verboseLevel>1) {
    G4cout << "hit in station "<<stationName<< " on technology "<<tech << G4endl;
    G4cout << "constructing ids (stName, stEta, stPhi, dr, dZ, dPhi)= "<<stationName<< " "<< stationEta<<" " << stationPhi<< " "<<doubletR<< " "<< doubletZ<< " "<<doubletPhi << G4endl;
  }
 
  HitID RPCid_eta = m_muonHelper->BuildRpcHitId(stationName, stationPhi, stationEta,
                                              mydbZ, doubletR, gasGap, mydbP,0);
 
  HitID RPCid_phi = m_muonHelper->BuildRpcHitId(stationName, stationPhi, stationEta,
                                              mydbZ, doubletR, gasGap, mydbP,1);
 
  // retrieve track barcode
  TrackHelper trHelp(aStep->GetTrack());
 
  //construct new rpc hit
  m_vertex = Amg::Hep3VectorToEigen( aStep->GetTrack()->GetVertexPosition() );
 
  // if the track m_vertex is far from (0,0,0), takes the tof, otherwise take the "usual" g4 globalTime
  (((m_vertex.mag()) < 100) ? (m_globalTime  = globalTime) : (m_globalTime = tof));
 
  m_myRPCHitColl->Emplace(RPCid_eta, m_globalTime,
                          localPosition,
                          trHelp.GenerateParticleLink(m_g4UserEventInfo ? m_g4UserEventInfo->GetEventStore() : nullptr),
                          localPostPosition,
                          aStep->GetTotalEnergyDeposit(),
                          aStep->GetStepLength(),
                          currentTrack->GetDefinition()->GetPDGEncoding(),
                          aStep->GetPreStepPoint()->GetKineticEnergy());
  m_myRPCHitColl->Emplace(RPCid_phi, m_globalTime,
                        localPosition,
                        trHelp.GenerateParticleLink(m_g4UserEventInfo ? m_g4UserEventInfo->GetEventStore() : nullptr),
                        localPostPosition,
                        aStep->GetTotalEnergyDeposit(),
                        aStep->GetStepLength(),
                        currentTrack->GetDefinition()->GetPDGEncoding(),
                        aStep->GetPreStepPoint()->GetKineticEnergy());
 
  return true;
}

Member Data Documentation

◆ m_currVertex

Amg::Vector3D RPCSensitiveDetectorCosmics::m_currVertex

private

Definition at line 121 of file RPCSensitiveDetectorCosmics.h.

◆ m_g4UserEventInfo

AtlasG4EventUserInfo* RPCSensitiveDetectorCosmics::m_g4UserEventInfo {nullptr}

private

Definition at line 114 of file RPCSensitiveDetectorCosmics.h.

114{nullptr};

◆ m_globalTime

double RPCSensitiveDetectorCosmics::m_globalTime

private

Definition at line 116 of file RPCSensitiveDetectorCosmics.h.

◆ m_globH

Amg::Vector3D RPCSensitiveDetectorCosmics::m_globH

private

Definition at line 122 of file RPCSensitiveDetectorCosmics.h.

◆ m_hitCollectionName

std::string RPCSensitiveDetectorCosmics::m_hitCollectionName

private

member data

Definition at line 112 of file RPCSensitiveDetectorCosmics.h.

◆ m_isGeoModel

bool RPCSensitiveDetectorCosmics::m_isGeoModel

private

Definition at line 117 of file RPCSensitiveDetectorCosmics.h.

◆ m_mom

Amg::Vector3D RPCSensitiveDetectorCosmics::m_mom

private

Definition at line 118 of file RPCSensitiveDetectorCosmics.h.

◆ m_momMag

double RPCSensitiveDetectorCosmics::m_momMag

private

Definition at line 119 of file RPCSensitiveDetectorCosmics.h.

◆ m_muonHelper

const RpcHitIdHelper* RPCSensitiveDetectorCosmics::m_muonHelper

private

Definition at line 115 of file RPCSensitiveDetectorCosmics.h.

◆ m_myRPCHitColl

RPCSimHitCollection* RPCSensitiveDetectorCosmics::m_myRPCHitColl {nullptr}

private

Definition at line 113 of file RPCSensitiveDetectorCosmics.h.

113{nullptr};

◆ m_vertex

Amg::Vector3D RPCSensitiveDetectorCosmics::m_vertex

private

Definition at line 120 of file RPCSensitiveDetectorCosmics.h.

The documentation for this class was generated from the following files:

Public Member Functions

Private Member Functions

Private Attributes

Detailed Description

Class methods and properties

notes:

Constructor & Destructor Documentation

◆ RPCSensitiveDetectorCosmics()

Member Function Documentation

◆ FRIEND_TEST() [1/2]

◆ FRIEND_TEST() [2/2]

◆ Initialize()

◆ ProcessHits()

Member Data Documentation

◆ m_currVertex

◆ m_g4UserEventInfo

◆ m_globalTime

◆ m_globH

◆ m_hitCollectionName

◆ m_isGeoModel

◆ m_mom

◆ m_momMag

◆ m_muonHelper

◆ m_myRPCHitColl

◆ m_vertex