g4simulation/g4main/PHG4EtaPhiParameterization.cc

0001 #include "PHG4EtaPhiParameterization.h"
0002
0003 #include <Geant4/G4ThreeVector.hh>
0004 #include <Geant4/G4Tubs.hh>
0005 #include <Geant4/G4Types.hh>  // for G4int
0006 #include <Geant4/G4VPhysicalVolume.hh>
0007
0008 #include <algorithm>  // for copy
0009 #include <cmath>
0010 #include <cstdlib>
0011 #include <iostream>
0012 #include <iterator>
0013
0014 PHG4EtaPhiParameterization::PHG4EtaPhiParameterization(
0015     unsigned int neta,  // Binning in eta
0016     double minEta,      // "
0017     double maxEta,      // "
0018     unsigned int nphi,  // number of phi bins
0019     double startPhi,
0020     double deltaPhi,
0021     double radiusIn,   // Radius of inner face of cylinder
0022     double radiusOut,  // Radius of outer face of cylinder
0023     double centerZ     // Z of center of set
0024     )
0025   : _neta(neta)
0026   , _minEta(minEta)
0027   , _maxEta(maxEta)
0028   , _nphi(nphi)
0029   , _startPhi(startPhi)
0030   , _deltaPhi(deltaPhi)
0031   , _radiusIn(radiusIn)
0032   , _radiusOut(radiusOut)
0033   , _centerZ(centerZ)
0034 {
0035   if (_maxEta < _minEta)
0036   {
0037     std::cout << " invalid eta, max<min"
0038               << " etamin: " << _minEta
0039               << " etamax: " << _maxEta
0040               << std::endl;
0041     exit(1);
0042   }
0043
0044   if ((_radiusIn < 0.0) || (_radiusOut < 0.0) || (_radiusOut < _radiusIn))
0045   {
0046     std::cout << " invalid radius parameters:"
0047               << " radiusIn: " << radiusIn
0048               << " radiusOut: " << radiusOut
0049               << std::endl;
0050     exit(1);
0051   }
0052
0053   double totalEta = _maxEta - _minEta;
0054   double dEta = totalEta / _neta;
0055   for (unsigned int i = 0; i < neta; i++)
0056   {
0057     // Compute the edges of this eta bin
0058     double etaMin = _minEta + dEta * i;
0059     double etaMax = etaMin + dEta;
0060     // double eta = _minEta + dEta * (i + 0.5);
0061     //  Compute the corresponding Z positions of the edges
0062     double zmin = _centerZ + _radiusIn * std::sinh(etaMin);
0063     double zmax = _centerZ + _radiusIn * std::sinh(etaMax);
0064     // Z positions is halfway between the edges
0065     double zpos = (zmin + zmax) / 2.0;
0066     double zhalf = (zmax - zmin) / 2.0;
0067     _zpos.push_back(zpos);
0068     _zhalf.push_back(zhalf);
0069   }
0070
0071   for (unsigned int i = 0; i < _nphi; i++)
0072   {
0073     _phi0.push_back(_startPhi + i * _deltaPhi);
0074     _phi1.push_back(_startPhi + (i + 1) * _deltaPhi);
0075   }
0076
0077   // Build lookup vectors for the copyNo->(ieta, iphi) translation
0078   //
0079   for (unsigned int i = 0; i < _neta * _nphi; i++)
0080   {
0081     div_t q = div((int) i, (int) _nphi);
0082     int ieta = q.quot;
0083     int iphi = q.rem;
0084     _ieta.push_back(ieta);
0085     _iphi.push_back(iphi);
0086   }
0087 }
0088
0089 PHG4EtaPhiParameterization::~PHG4EtaPhiParameterization()
0090 {
0091   std::cout << "PHG4EtaPhiParameterization::~PHG4EtaPhiParameterization: Alas, poor Yorick! I knew him, Horatio"
0092             << std::endl;
0093 }
0094
0095 void PHG4EtaPhiParameterization::Print(std::ostream& os) const
0096 {
0097   os << "PhiEtaPhiParameterization: NETA x NPHI = " << _neta << " x " << _nphi << std::endl;
0098   os << "Zpos: ";
0099   std::copy(_zpos.begin(), _zpos.end(), std::ostream_iterator<double>(os, " "));
0100   os << std::endl;
0101   os << "Phi0: ";
0102   std::copy(_phi0.begin(), _phi0.end(), std::ostream_iterator<double>(os, " "));
0103   os << std::endl;
0104   os << "Phi1: ";
0105   std::copy(_phi0.begin(), _phi0.end(), std::ostream_iterator<double>(os, " "));
0106   os << std::endl;
0107 }
0108
0109 void PHG4EtaPhiParameterization::ComputeTransformation(const G4int copyNo, G4VPhysicalVolume* physVol) const
0110 {
0111   int iring = copyNo / _nphi;
0112   G4ThreeVector origin(0, 0, _zpos.at(iring));
0113   physVol->SetTranslation(origin);
0114   physVol->SetRotation(nullptr);
0115 }
0116
0117 void PHG4EtaPhiParameterization::ComputeDimensions(G4Tubs& ring, const G4int copyNo,
0118                                                    const G4VPhysicalVolume* /*unused*/) const
0119 {
0120   int ieta = GetIEta(copyNo);
0121   int iphi = GetIPhi(copyNo);
0122   double phi = _phi0.at(iphi);
0123   double dphi = _phi1.at(iphi) - phi;
0124   ring.SetZHalfLength(_zhalf.at(ieta));
0125   ring.SetInnerRadius(_radiusIn);
0126   ring.SetOuterRadius(_radiusOut);
0127   ring.SetStartPhiAngle(phi);
0128   ring.SetDeltaPhiAngle(dphi);
0129 }