DiffusionGTS.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/DiffusionGTS.h"
 
#include "TRandom.h"
 
#include <iostream>
#include <fstream>
#include <math.h>
 
using namespace std;
 
DiffusionGTS::DiffusionGTS(ICgemGeomSvc* geomSvc) {
  m_geomSvc = geomSvc;
}
DiffusionGTS::~DiffusionGTS(){}
 
void DiffusionGTS::readGasPerpParameters(std::string fileName){
 
  ifstream fin(fileName.c_str(), ios::in);
  
  string strline;
  string strpar;
  if( ! fin.is_open() ){
    std::cout << "DiffusionGTS::readGasPerpParameters ERROR: can not open file " << fileName.c_str() << endl;
  }
  while(fin >> strpar){
    if(strpar[0] == '#'){
      getline(fin, strline);
    }
    else if(strpar == "n_shift_x_drift"){
      fin >> n_shift_x_drift;
    }
    else if(strpar == "shift_x_drift_par"){
      for(int i=0; i<n_shift_x_drift; i++) fin >> m_shift_x_drift_par[i]; 
    }
    else if(strpar == "n_sigma_x_drift"){
      fin >> n_sigma_x_drift;
    }
    else if(strpar == "sigma_x_drift_par"){
      for(int i=0; i<n_sigma_x_drift; i++) fin >> m_sigma_x_drift_par[i];
    }
    else if(strpar == "shift_x_transf_par"){
      fin >> m_shift_x_transf_par;
    }
    else if(strpar == "sigma_x_transf_par"){
      fin >> m_sigma_x_transf_par;
    }
    else if(strpar == "shift_x_induct_par"){
      fin >> m_shift_x_induct_par;
    }
    else if(strpar == "sigma_x_induct_par"){
      fin >> m_sigma_x_induct_par;
    }
    else if(strpar == "n_shift_t_drift"){
      fin >> n_shift_t_drift;
    }
    else if(strpar == "shift_t_drift_par"){
      for(int i=0; i<n_shift_t_drift; i++) fin >> m_shift_t_drift_par[i];
    }
    else if(strpar == "n_sigma_t_drift"){
      fin >> n_sigma_t_drift;
    }
    else if(strpar == "sigma_t_drift_par"){
      for(int i=0; i<n_sigma_t_drift; i++) fin >> m_sigma_t_drift_par[i];
    }
    else if(strpar == "shift_t_transf_par"){
      fin >> m_shift_t_transf_par;
    }
    else if(strpar == "sigma_t_transf_par"){
      fin >> m_sigma_t_transf_par;
    }
    else if(strpar == "shift_t_induct_par"){
      fin >> m_shift_t_induct_par;
    }
    else if(strpar == "sigma_t_induct_par"){
      fin >> m_sigma_t_induct_par;
    }
  }
 
  std::cout << "DiffusionGTS::readGasPerpParameters: open file " << fileName.c_str() << endl;
 
  //  cout << n_shift_x_drift  << endl;
  cout << "shift x in drift";  
  for(int i=0; i<n_shift_x_drift; i++) cout << " " << m_shift_x_drift_par[i];
  cout << endl;
  // cout << n_sigma_x_drift << endl;
  cout << "sigma x in drift";
  for(int i=0; i<n_sigma_x_drift; i++) cout << " " << m_sigma_x_drift_par[i];
  cout << endl;
  cout << "shift x in transf " << m_shift_x_transf_par << endl;
  cout << "sigma x in transf " << m_sigma_x_transf_par << endl;
  cout << "shift x in induct " << m_shift_x_induct_par << endl;
  cout << "sigma x in induct " << m_sigma_x_induct_par << endl;
 
  cout << "shift t in drift";
  // cout << n_shift_t_drift << endl;
  for(int i=0; i<n_shift_t_drift; i++) cout << " " << m_shift_t_drift_par[i];
  cout << endl;
  //  cout << n_sigma_t_drift << endl;
  cout << "sigma t in drift";
  for(int i=0; i<n_sigma_t_drift; i++) cout << " " << m_sigma_t_drift_par[i];
  cout << endl;
  cout << "shift t in transf " << m_shift_t_transf_par << endl;
  cout << "sigma t in transf " << m_sigma_t_transf_par << endl;
  cout << "shift t in induct " << m_shift_t_induct_par << endl;
  cout << "sigma t in induct " << m_sigma_t_induct_par << endl;
  
}
// --- WITHOUT FIELD ---
void DiffusionGTS::readGasParalParameters(std::string fileName){
 
  ifstream fin(fileName.c_str(), ios::in);
  
  string strline;
  string strpar;
  if( ! fin.is_open() ){
    std::cout << "DiffusionGTS::readGasParalParameters ERROR: can not open file " << fileName.c_str() << endl;
  }
  while(fin >> strpar){
    if(strpar[0] == '#'){
      getline(fin, strline);
    }
    else if(strpar == "n_shift_x_drift"){
      fin >> n_shift_y_drift;
    }
    else if(strpar == "shift_x_drift_par"){
      for(int i=0; i<n_shift_y_drift; i++) fin >> m_shift_y_drift_par[i]; 
    }
    else if(strpar == "n_sigma_x_drift"){
      fin >> n_sigma_y_drift;
    }
    else if(strpar == "sigma_x_drift_par"){
      for(int i=0; i<n_sigma_y_drift; i++) fin >> m_sigma_y_drift_par[i];
    }
    else if(strpar == "shift_x_transf_par"){
      fin >> m_shift_y_transf_par;
    }
    else if(strpar == "sigma_x_transf_par"){
      fin >> m_sigma_y_transf_par;
    }
    else if(strpar == "shift_x_induct_par"){
      fin >> m_shift_y_induct_par;
    }
    else if(strpar == "sigma_x_induct_par"){
      fin >> m_sigma_y_induct_par;
    }
    else break;
  }
 
 
  std::cout << "DiffusionGTS::readGasParalParameters: open file " << fileName.c_str() << endl;
 
  cout << "shift y in drift";  
  for(int i=0; i<n_shift_y_drift; i++) cout << " " << m_shift_y_drift_par[i];
  cout << endl;
  cout << "sigma y in drift";
  for(int i=0; i<n_sigma_y_drift; i++) cout << " " << m_sigma_y_drift_par[i];
  cout << endl;
  cout << "shift y in transf " <<  m_shift_y_transf_par << endl;
  cout << "sigma y in transf " << m_sigma_y_transf_par << endl;
  cout << "shift y in induct " << m_shift_y_induct_par << endl;
  cout << "sigma y in induct " << m_sigma_y_induct_par << endl;
 
  cout << "shift t in drift";
  // cout << n_shift_t_drift << endl;
  for(int i=0; i<n_shift_t_drift; i++) cout << " " << m_shift_t_drift_par[i];
  cout << endl;
  
}
// ----------------- X -------------------
double DiffusionGTS::shift_x_drift(double z){
 
  double dx = 0;  
  /**
     if(z<0||z>5.) {  // CHECK 16) gap from geometry! 
     cout<<"DiffusionGTS::shift_x_drift: error in z = " << z << endl;
     return dx;
     } 
  **/
  for(int i=0; i<n_shift_x_drift; i++) dx += pow(z, i) * m_shift_x_drift_par[i];
  return dx;
}
double DiffusionGTS::shift_x_transf() {
  double dx = m_shift_x_transf_par;
  return dx;
}
double DiffusionGTS::shift_x_induct() {
  double dx = m_shift_x_induct_par;
  return dx;
}
 
 
double DiffusionGTS::sigma_x_drift(double z){
  double dx = 0;  
  /**
     if(z<0||z>5.) {  // CHECK 16) gap from geometry! 
     cout<<"DiffusionGTS::sigma_x_drift: error in z = " << z << endl;
     return dx;
     } 
  **/
  for(int i=0; i<n_sigma_x_drift; i++) dx += pow(z, i) * m_sigma_x_drift_par[i];
  return dx;
}
double DiffusionGTS::sigma_x_transf() {
  double dx = m_sigma_x_transf_par;
  return dx;
}
double DiffusionGTS::sigma_x_induct() {
  double dx = m_sigma_x_induct_par;
  return dx;
}
 
//----------------- Y
double DiffusionGTS::shift_y_drift(double z){
 
  double dy = 0;  
  /**
     if(z<0||z>5.) {  // CHECK 16) gap from geometry! 
     cout<<"DiffusionGTS::shift_y_drift: error in z = " << z << endl;
     return dy;
     } 
  **/
  for(int i=0; i<n_shift_y_drift; i++) dy += pow(z, i) * m_shift_y_drift_par[i];
  return dy;
}
double DiffusionGTS::shift_y_transf() {
  double dy = m_shift_y_transf_par;
  return dy;
}
double DiffusionGTS::shift_y_induct() {
  double dy = m_shift_y_induct_par;
  return dy;
}
 
double DiffusionGTS::sigma_y_drift(double z){
  double dy = 0;  
  /**
     if(z<0||z>5.) {  // CHECK 16) gap from geometry! 
     cout<<"DiffusionGTS::sigma_y_drift: error in z = " << z << endl;
     return dy;
     } 
  **/
  for(int i=0; i<n_sigma_y_drift; i++) dy += pow(z, i) * m_sigma_y_drift_par[i];
  return dy;
}
double DiffusionGTS::sigma_y_transf() {
  double dy = m_sigma_y_transf_par;
  return dy;
}
double DiffusionGTS::sigma_y_induct() {
  double dy = m_sigma_y_induct_par;
  return dy;
}
 
// ----------------------------------------
double DiffusionGTS::shift_t_drift(double z){
  double dt = 0;
  /**
     if(z<0||z>5.) {  // CHECK 16) gap from geometry!
     cout<<"DiffusionGTS::shift_t_drift: error in z = " << z << endl;
     return dt;
     }
  **/
  for(int i=0; i<n_shift_t_drift; i++) dt += pow(z, i) * m_shift_t_drift_par[i];
  return dt;
}
 
double DiffusionGTS::shift_t_transf() {
  double dx = m_shift_t_transf_par;
  return dx;
}
double DiffusionGTS::shift_t_induct() {
  double dx = m_shift_t_induct_par;
  return dx;
}
 
double DiffusionGTS::sigma_t_drift(double z){
  double dt = 0;
  /**
     if(z<0||z>5.) {  // CHECK 16) gap from geometry! 
     cout<<"DiffusionGTS::sigma_t_drift: error in z = " << z << endl;
     return dt;
     }
  **/
  for(int i=0; i<n_sigma_t_drift; i++) dt += pow(z, i) * m_sigma_t_drift_par[i];
  return dt;
}
double DiffusionGTS::sigma_t_transf() {
  double dx = m_sigma_t_transf_par;
  return dx;
}
double DiffusionGTS::sigma_t_induct() {
  double dx = m_sigma_t_induct_par;
  return dx;
}
 
 
 
// ----------------------------------------------------------- drifting
void DiffusionGTS::diffusion_drift(double xi, double yi, double zi, double ti, double &xf, double &yf, double &zf, double &tf) {
 
  // x
  double dx = shift_x_drift(zi);
  double sx = sigma_x_drift(zi);
  xf = xi + gRandom->Gaus(dx, sx);
 
  // y
  double dy = shift_y_drift(zi);
  double sy = sigma_y_drift(zi);
  yf = yi + gRandom->Gaus(dy, sy);
 
  // t
  double dt = shift_t_drift(zi);
  double st = sigma_t_drift(zi);
  tf = ti + gRandom->Gaus(dt, st);
}
 
void DiffusionGTS::diffusion_transf(double xi, double yi, double zi, double ti, double &xf, double &yf, double &zf, double &tf) {
 
  // x
  double dx = shift_x_transf();
  double sx = sigma_x_transf();
  xf = xi + gRandom->Gaus(dx, sx);
 
  // y
  double dy = shift_y_transf();
  double sy = sigma_y_transf();
  yf = yi + gRandom->Gaus(dy, sy);
 
  // t
  double dt = shift_t_transf();
  double st = sigma_t_transf();
  tf = ti + gRandom->Gaus(dt, st);
}
 
 
void DiffusionGTS::diffusion_induct(double xi, double yi, double zi, double ti, double &xf, double &yf, double &zf, double &tf) {
 
  // x
  double dx = shift_x_induct();
  double sx = sigma_x_induct();
  xf = xi + gRandom->Gaus(dx, sx);
 
  // y
  double dy = shift_y_induct();
  double sy = sigma_y_induct();
  yf = yi + gRandom->Gaus(dy, sy);
 
  // t                      
  double dt = shift_t_induct();
  double st = sigma_t_induct();
  tf = ti + gRandom->Gaus(dt, st);
 
}
 
// -------------------------------------