APV25GTS.cxx

Go to the documentation of this file.

// implementation: 
//          R. Farinelli (University of Ferrara & INFN of Ferrara)
//          L. Lavezzi   (IHEP & INFN of Torino)
#include "CgemDigitizerSvc/APV25GTS.h"
 
#include "TRandom.h"
#include "TStyle.h"
#include "TF1.h"
#include "TMath.h"
 
#define ele_to_fC 1.6e-4 // CHECK
#define tstep 25         // CHECK
 
APV25GTS::APV25GTS() : histo_q_vs_t(NULL) {
  histo_q_vs_t = new TH1F("histo_q_vs_t", "", 27, 0, 675); // 25 ns x 27 time bins
}
 
APV25GTS::~APV25GTS(){
  delete histo_q_vs_t;
}
 
void APV25GTS::add_charge(double q, double t) {
  histo_q_vs_t->Fill(t, q);
}
 
 
double APV25GTS::compute_charge() {              // in fC
  double qmax = histo_q_vs_t->GetMaximum();
  double qmax_fC = qmax/ele_to_fC;
  return qmax_fC;
}
 
double APV25GTS::compute_tFD() {
 
  // TIME from FERMI DIRAC %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%                                                                        
  // start parameters                                                                                                                     
  double maxbin    = histo_q_vs_t->GetMaximumBin();
  double QMaxHisto = compute_charge();
 
  // double startBin  = maxbin - 5;                                                                                                       
  // double endBin    = maxbin + 10;                                                                                                      
  //find the lower edge of the rising @ 10% of the maximum                                                                                
  double minbin = 0;
  for(int ibin = maxbin; ibin > 1; ibin--){
    if(histo_q_vs_t->GetBinContent(ibin) < 0.1 * QMaxHisto) {
      minbin = ibin;
      break;
    }
  }
  double startvalues[5] = {0};
  double fitlimup[5]    = {0};
  double fitlimlow[5]   = {0};
  // media dei tre bin prima del minbin                                                                                                   
  double meanfirstbin = (histo_q_vs_t->GetBinContent(minbin-1) + histo_q_vs_t->GetBinContent(minbin-2) + histo_q_vs_t->GetBinContent(minbin-3))/3.;
  double sigma_MFB = TMath::Sqrt( ( pow(histo_q_vs_t->GetBinContent(minbin-1) - meanfirstbin, 2) + pow(histo_q_vs_t->GetBinContent(minbin-2) - meanfirstbin,2) + pow(histo_q_vs_t->GetBinContent(minbin-3)- meanfirstbin,2)) / 3.);
  if(sigma_MFB < TMath::Abs(meanfirstbin)) sigma_MFB = TMath::Abs(meanfirstbin);
 
  startvalues[0] = meanfirstbin;
  startvalues[1] = QMaxHisto;
  startvalues[2] = 0.5* (minbin + maxbin) * tstep;
  startvalues[3] = 0.5*tstep;
 
  // put liimits                                                                                                                          
  fitlimlow[0] = startvalues[0] - 2 * sigma_MFB;
  fitlimlow[1] = startvalues[1] - 0.3 * startvalues[1];
  fitlimlow[2] = minbin * tstep;
  fitlimlow[3] = 0.1*tstep;
 
  fitlimup[0] = startvalues[0] + 2 * sigma_MFB;
  fitlimup[1] = startvalues[1] + 0.3 * startvalues[1];
  fitlimup[2] = maxbin * tstep;
  fitlimup[3] = 0.9*tstep;
  // ...................                                                                                                                  
  // the real Fermi Dirac                                                                                                                 
  //                        1                                                                                                             
  // FD(t) = K ----------------------------- + B                                                                                          
  //            1 + exp( -(t - tFD)/sigFD )                                                                                               
  // [0] = B                                                                                                                              
  // [1] = K                                                                                                                              
  // [2] = tFD                                                                                                                            
  // [3] = sigFD                                                                                                                          
  double minFD = minbin - 4;
  double maxFD = maxbin + 0;
 
  TF1 f_FD("f_FD", "[0] + [1]/(1+TMath::Exp(-(x - [2])/[3]))", minFD*tstep, maxFD*tstep);
  f_FD.SetParameters(startvalues[0], startvalues[1], startvalues[2], startvalues[3]);
  f_FD.SetParLimits(0, fitlimlow[0], fitlimup[0]);
  f_FD.SetParLimits(1, fitlimlow[1], fitlimup[1]);
  f_FD.SetParLimits(2, fitlimlow[2], fitlimup[2]);
  f_FD.SetParLimits(3, fitlimlow[3], fitlimup[3]);
 
  gStyle->SetOptFit(1111);
  histo_q_vs_t->Fit("f_FD","WQRB");  // quiet/range/params                                                                            
 
  double tFD = 0;
  tFD = gRandom->Gaus(f_FD.GetParameter(2), 8);
  double dtFD = 0.;
  dtFD = f_FD.GetParError(2);
  set_dtFD(dtFD);
 
  // save parameters                                                                                                                      
  set_FD_par(0, f_FD.GetParameter(0));
  set_FD_par(1, f_FD.GetParameter(1));
  set_FD_par(2, f_FD.GetParameter(2));
  set_FD_par(3, f_FD.GetParameter(3));
 
 
  // CHECK
  /**
  double zzz = (tFD-strip_t0-Get_Track_T0())*drift_velocity; // chech --> // This Z assign Z=0.5 to the cathode and Z=0 to G1 CHECK      
  set_APV25_Z(zzz);
 
  if(Get_APV25_Q(iS)>qcut && Get_APV25_T(iS)!=0){
    Set_APV25_X(iS,iS*pitch+0.5*pitch+x_min);
  }
  **/
 
  return tFD;
}