CbmTrdParModGas.cxx

Go to the documentation of this file.

#include "CbmTrdParModGas.h"
 
#include <FairLogger.h>  // for Logger, LOG
 
#include <TAxis.h>              // for TAxis
#include <TGenericClassInfo.h>  // for TGenericClassInfo
#include <TH2.h>                // for TH2F
#include <TMath.h>              // for Exp
#include <TMathBase.h>          // for Abs
#include <TObjArray.h>          // for TObjArray
#include <TObjString.h>         // for TObjString
 
#include <stdio.h>  // for printf
 
#define VERBOSE 0
 
// binding energy in keV for 'K' 'L' and 'M' shells of Ar and Xe
Float_t CbmTrdParModGas::fgkBindingEnergy[2][NSHELLS] = {
  {34.5, 5.1, 1.1},  // Xe
  {3.2, 0.25, 0.1}   // Ar
};
// binding energy in keV for 'K' 'L' and 'M' shells of Ar and Xe
Float_t CbmTrdParModGas::fgkBR[2][NSHELLS - 1] = {
  {0.11, 0.02},  // Xe
  {0.11, 0.02}   // Ar
};
// average energy to produce one electron-ion pair
Float_t CbmTrdParModGas::fgkWi[3] = {
  23.1,  // Xe
  25.8,  // Ar
  33     // CO2
};
// Bucharest detector gas gain parametrization based on 55Fe measurements with ArCO2(80/20)
Float_t CbmTrdParModGas::fgkGGainUaPar[2] = {-10.1676, 8.3745};
Float_t CbmTrdParModGas::fgkE0 = 866.1047;  // energy offset in ADC ch @ 0 keV
 
//___________________________________________________________________
CbmTrdParModGas::CbmTrdParModGas(const char* title)
  : CbmTrdParMod("CbmTrdParModGas", title)
  , fConfig(0)
  , fUa(0)
  , fUd(0)
  , fDw(0.3)
  , fGasThick(0.6)
  , fPercentCO2(0.2)
  , fDriftMap(nullptr)
  , fFileNamePID() {
  TString s(title);
  TString name;
  //  Int_t val;
  TObjArray* so = s.Tokenize("/");
  for (Int_t ie(0); ie < so->GetEntries(); ie += 2) {
    name = ((TObjString*) (*so)[ie])->String();
    if (name.EqualTo("Module"))
      fModuleId = ((TObjString*) (*so)[ie + 1])->String().Atoi();
    else if (name.EqualTo("Ua"))
      fUa = ((TObjString*) (*so)[ie + 1])->String().Atoi();
    else if (name.EqualTo("Ud"))
      fUd = ((TObjString*) (*so)[ie + 1])->String().Atoi();
    else if (name.EqualTo("Gas")) {
      TString gas = ((TObjString*) (*so)[ie + 1])->String();
      if (gas.EqualTo("Ar"))
        SetNobleGasType(1);
      else if (gas.EqualTo("Xe"))
        SetNobleGasType(0);
      else {
        LOG(warn) << GetName() << ":: gas type \"" << gas
                  << "\" not defined. Default to Xe.";
        SetNobleGasType(0);
      }
    }
  }
  so->Delete();
  delete so;
 
  if (VERBOSE) printf("Module[%2d] U[%4d %3d]\n", fModuleId, fUa, fUd);
}
 
//___________________________________________________________________
CbmTrdParModGas::~CbmTrdParModGas() {
  //  if(fDriftMap) delete fDriftMap;
}
 
//_______________________________________________________________________________________________
Float_t CbmTrdParModGas::EkevFC(Float_t ekev) const {
  // Convert energy deposit to no of primary ionizations and apply gas gain.
  // Currently gas gain is evalauted from 55Fe spectrum analysis on ArCO2(80/20)
  Int_t gasId = GetNobleGasType() - 1;
  Float_t wi  = (1. - fPercentCO2) * fgkWi[gasId] + fPercentCO2 * fgkWi[2];
 
  //gas gain
  // G = G[ev->ADC] * wi[ArCO2]/C[mV->ADC]/A[fC->mV]/e
  // G[ev->ADC] : measured gain based on 55Fe spectrum (expo)
  // wi[ArCO2]  : average energy to produce a ele-ion pair in mixture (27.24 ev)
  // C[mV->ADC] : FASP out [2V] to ADC range [4096 ch] (2)
  // A[fC->mV]  : FASP gain from CADENCE (6)
  // e : 1.6e-4 [fC] electric charge
  Double_t gain =
    170.25 * TMath::Exp(fgkGGainUaPar[0] + fgkGGainUaPar[1] * fUa * 1.e-3)
    / 12.;
  // for Xe correct Ar gain measurements TODO
  if (gasId == 0) gain *= 0.6;
 
  Float_t efC = gain * ekev * 0.16 / wi;
  if (VERBOSE)
    printf(
      "        ua[V]=%d gain[%5.2e] wi[eV]=%5.2f :: E[keV]=%6.3f E[fC]=%6.2f\n",
      fUa,
      gain,
      wi,
      ekev,
      efC);
  return efC;
 
  //   Double_t eadc = fgkE0 + gain * ekev;
  //   //printf("E = %fkeV %fADC\n", ekev, eadc);
  //   // apply FASP gain -> should be done without intermediate ADC conversion TODO
  //   Double_t sFASP = eadc/2.; // FASP signal [mV]; FASP amplification 1
  //   // FASP gaincharacteristic -> should be defined elsewhere
  //   // data based on CADENCE simulations
  //   Double_t s0FASP = 10, gFASP = 6;
  //   return (sFASP-s0FASP)/gFASP;
}
 
//_______________________________________________________________________________________________
Int_t CbmTrdParModGas::GetShellId(const Char_t shell) const {
  switch (shell) {
    case 'K': return 0;
    case 'L': return 1;
    case 'M': return 2;
    default:
      LOG(warn) << GetName() << "::GetShellId: Atomic shell : " << shell
                << " not defined for gas "
                << (GetNobleGasType() == 2 ? "Ar" : "Xe");
      return -1;
  }
}
 
//_______________________________________________________________________________________________
Float_t CbmTrdParModGas::GetBindingEnergy(const Char_t shell,
                                          Bool_t main) const {
  Int_t gasId   = GetNobleGasType() - 1;
  Int_t shellId = GetShellId(shell);
  if (shellId < 0) return 0;
 
  if (!main)
    return fgkBindingEnergy[gasId][shellId];
  else {
    if (shellId < NSHELLS - 1)
      return fgkBindingEnergy[gasId][shellId + 1];
    else {
      LOG(warn) << GetName()
                << "::GetBindingEnergy: Request atomic shell : " << shellId + 1
                << " not defined for gas " << (gasId ? "Ar" : "Xe");
      return 0;
    }
  }
  return 0;
}
 
//_______________________________________________________________________________________________
Float_t CbmTrdParModGas::GetNonIonizingBR(const Char_t shell) const {
  Int_t gasId   = GetNobleGasType() - 1;
  Int_t shellId = GetShellId(shell);
  if (shellId < 0) return 0;
 
  return fgkBR[gasId][shellId];
}
 
//_______________________________________________________________________________________________
Char_t CbmTrdParModGas::GetPEshell(Float_t Ex) const {
  const Char_t shellName[NSHELLS] = {'K', 'L', 'M'};
  Int_t gasId                     = GetNobleGasType() - 1;
  for (Int_t ishell(0); ishell < NSHELLS; ishell++) {
    if (Ex < fgkBindingEnergy[gasId][ishell]) continue;
    return shellName[ishell];
  }
  LOG(debug) << GetName() << "::GetPEshell: Ex[keV] " << Ex
             << " less than highes atomic shell binding energy : "
             << fgkBindingEnergy[gasId][NSHELLS - 1] << " for gas "
             << (gasId ? "Ar" : "Xe");
  return 0;
}
 
//_______________________________________________________________________________________________
Double_t CbmTrdParModGas::GetDriftTime(Double_t y0, Double_t z0) const {
  const TAxis *ay(fDriftMap->GetXaxis()), *az(fDriftMap->GetYaxis());
  Int_t by(ay->FindBin(y0)), bz(az->FindBin(z0));
  Double_t tmin(fDriftMap->GetBinContent(by, bz));
  if (VERBOSE) printf("GetDriftTime :: Start @ dt=%3d [ns]\n", Int_t(tmin));
  return tmin;
}
 
//_______________________________________________________________________________________________
void CbmTrdParModGas::Print(Option_t* /*opt*/) const {
  printf("%s @ %4d ", GetName(), fModuleId);
  printf("Type[%s] ", GetDetName());
  printf("%s[%4.1f%%] Ua[V]=%d Ud[V]=%d ",
         GetNobleGasName(),
         1e2 * GetNobleGas(),
         fUa,
         fUd);
  printf("Pid Type[%d] DB[%s]\n", GetPidType(), fFileNamePID.Data());
}
 
//_______________________________________________________________________________________________
Double_t CbmTrdParModGas::ScanDriftTime(Double_t y0,
                                        Double_t z0,
                                        Double_t dzdy,
                                        Double_t dy) const {
  Double_t y1 = y0 + dy, z1 = z0 + dzdy * dy, dw(fDw), dwh(0.5 * dw);
  //  Double_t dhh(fGasThick);
 
  if (VERBOSE)
    printf("ScanDriftTime :: Try : [%7.4f %7.4f] => [%7.4f %7.4f] dzdy[%5.2f] "
           "dy[%5.2f]\n",
           y0,
           z0,
           y1,
           z1,
           dzdy,
           dy);
  while (y1 < -dwh - 1e-3) {
    y0 += dw;
    y1 += dw;
  }
  while (y1 > dwh + 1.e-3) {
    y0 -= dw;
    y1 -= dw;
  }
  TH2F* h = fDriftMap;
  y1      = y0;
  z1      = z0;
  const TAxis *ay(fDriftMap->GetXaxis()), *az(fDriftMap->GetYaxis());
  Int_t by(ay->FindBin(y1)), bz(az->FindBin(z1)), nby(ay->GetNbins()),
    nbz(az->GetNbins());
  Float_t dyStep = ay->GetBinWidth(1), tmin(500), tmax(0), tc(0);
  while (by > 0 && by <= nby && bz <= nbz) {
    bz = az->FindBin(z1);
    tc = h->GetBinContent(by, bz);
    //if(VERBOSE) printf("ScanDriftTime :: Do  : y0(%7.4f), z0(%7.4f), by(%3d), bz(%3d) td[ns]=%5.1f  t0[ns]=%5.1f \n", y1, z1, by, bz, tc, tmin);
 
    if (tc > 0) {
      if (tc < tmin) tmin = tc;
      if (tc > tmax) tmax = tc;
    }
    z1 += TMath::Abs(dzdy) * dyStep;
    if (dzdy > 0) {
      y1 += dyStep;
      by++;
    } else {
      y1 -= dyStep;
      by--;
    }
  }
  if (VERBOSE) printf("ScanDriftTime :: Start @ dt=%3d [ns]\n", Int_t(tmin));
 
  return tmin;
}
 
//___________________________________________________________________
void CbmTrdParModGas::SetDriftMap(TH2F* hm, TDirectory* d) {
  if (VERBOSE)
    printf("CbmTrdParModGas::SetDriftMap : Module[%2d] U[%4d %3d]\n",
           fModuleId,
           fUa,
           fUd);
 
  fDriftMap = (TH2F*) hm->Clone(Form("trdDM%02d", fModuleId));
  fDriftMap->SetTitle(GetTitle());
  fDriftMap->SetDirectory(d);
}
 
ClassImp(CbmTrdParModGas)