#include "PamVMCTrkDig.h"
#include <TMath.h>

ClassImp(PamVMCTrkDig)

  using TMath::Abs;

void PamVMCTrkDig::LoadFile(){

  cout<<"Loading Tracker Calibrations..."<<endl; 

  Int_t time=0;
  Int_t type=8;

  fdberr = fsql->Query_GL_PARAM(time,type); 
  //we hardcoded values of time and type of calibrations
  //later we'll introduce an options which allow us to load
  //from different sources
  
  if(fdberr<0){
    TString crfname = fpath+"/CalibTrk_00110_000_000.root"; 
    cout<<"No such record in DB for TRK: time="<<time
	<<" type="<<type<<endl
	<<"We will use file:"<<crfname<<endl;

    fcrfile = new TFile(crfname);

  } else {

    fquery.str("");
    fquery << fsql->GetPAR()->PATH.Data() << "/";
    fquery << fsql->GetPAR()->NAME.Data();

    ThrowCalFileUsage("TRK",fquery.str().c_str());

    fcrfile = new TFile(fquery.str().c_str());
  }

  if(!fcrfile) ThrowCalFileWarning("TRK"); else
    if(fcrfile->IsZombie()){
      ThrowCalFileWarning("TRK"); //critical error
      return;
    }

  //Load calibrations. I don't know, some empiric numbers of
  //broken VA chips introduced. Should be load as options..

            //------X-plane------//
  TTree *tr1 = (TTree*)fcrfile->Get("CalibTrk1");
  if ( !tr1 ){
    cout<<"!!!WARNING Tree CalibTrk1 in file "<<fcrfile->GetName()
        <<" was NOT found!!!"<<endl;
    return;
  }
  CalibTrk1Event *caldata1 = 0;
  TBranch *trbr1 = tr1->GetBranch("CalibTrk1");
  tr1->SetBranchAddress("CalibTrk1", &caldata1);
  trbr1->GetEntry(0);

  Int_t Kview,jj;

  for (Int_t i=0; i<6;i++) {
    for (Int_t j=0; j<3072;j++) {
      jj=j;
      // broken va1 replaced
      BrokenStrip(0,4,i,j,jj);
      BrokenStrip(5,12,i,j,jj);

      Kview=caldata1->DSPnumber[i]-1;
      fPedeTrack[Kview][j]=caldata1->DSPped_par[i][jj];
      fSigmaTrack[Kview][j]=caldata1->DSPsig_par[i][jj];
      if(caldata1->DSPbad_par[i][jj]==1) fSigmaTrack[Kview][j]=-fSigmaTrack[Kview][j];
    };
  };

            //------Y-plane------//
  TTree *tr2 = (TTree*)fcrfile->Get("CalibTrk2");
  if ( !tr2 ){
    cout<<"!!!WARNING Tree CalibTrk2 in file "<<fcrfile->GetName()
        <<" was NOT found!!!"<<endl;
    return;
  }
  CalibTrk2Event *caldata2 = 0;
  TBranch *trbr2 = tr2->GetBranch("CalibTrk2");
  tr2->SetBranchAddress("CalibTrk2",&caldata2);
  trbr2->GetEntry(1);

  for (Int_t i=0; i<6;i++) {
    for (Int_t j=0; j<3072;j++) {
      jj=j;
      // broken va1 replaced
      BrokenStrip(2,3,i,j,jj);
      BrokenStrip(2,5,i,j,jj);
      BrokenStrip(2,6,i,j,jj);
      BrokenStrip(2,7,i,j,jj);
      BrokenStrip(2,11,i,j,jj);

      Kview=caldata2->DSPnumber[i]-1;
      fPedeTrack[Kview][j]=caldata2->DSPped_par[i][jj];
      fSigmaTrack[Kview][j]=caldata2->DSPsig_par[i][jj];
      if(caldata2->DSPbad_par[i][jj]==1) fSigmaTrack[Kview][j]=-fSigmaTrack[Kview][j];  
    };
  };

  fcrfile->Close();
}


void PamVMCTrkDig::DigitizeTrackCalib(Int_t n){
  
  cout<<"Starting Tracker Calibrations..."<<endl;
  if( (n!=1)&&(n!=2) ) {
    cout << "!!!ERROR: Wrong DigitizeTrackCalib argument!!!" << endl;
    return;
  }; 

  fData.clear();

  UShort_t Dato;

  USBuffer tempdat; //temporary buffer to keep data before writing DSP data
  USBuffer::const_iterator p; //iterator;

  Float_t  dato1,dato2,dato3,dato4;
  UShort_t DatoDec,DatoDec1,DatoDec2,DatoDec3,DatoDec4;
  UShort_t EVENT_CAL, PED_L1, ReLength,OveCheckCode;
  UShort_t CkSum;

  for (Int_t j=n-1; j<fNviews;j+=2) { //0(1)...12
    CkSum=0;
    // here skip the dsp header and his trailer , to be written later
    tempdat.clear(); //clearing temporary buffer
    for (Int_t i=0; i<fNladder;i++) { //0...2
      for (Int_t k=0; k<fNstrips_ladder;k++) {  //0...1023
	// write in buffer the current LADDER
	Dato=(UShort_t)fPedeTrack[j][i*fNstrips_ladder+k];
	dato1=fPedeTrack[j][i*fNstrips_ladder+k]-Dato; 

	DatoDec1=(UShort_t)(dato1*2);
	dato2=dato1*2-DatoDec1;

	DatoDec2=(UShort_t)(dato2*2);
	dato3=dato2*2-DatoDec2;
      
	DatoDec3=(UShort_t)(dato3*2);
	dato4=dato3*2-DatoDec3;
      
	DatoDec4=(UShort_t)(dato4*2);
      
	DatoDec=DatoDec1*0x0008+DatoDec2*0x0004+DatoDec3*0x0002+DatoDec4*0x0001;

	tempdat.push_back( ((Dato << 4) | (DatoDec & 0x000F)) );
	CkSum^=tempdat.back();
      };

      for (Int_t k=0; k<fNstrips_ladder;k++) { //0...1023
      // write in buffer the current LADDER
	Dato=(UShort_t)fabs(fSigmaTrack[j][i*fNstrips_ladder+k]); 
	dato1=fabs(fSigmaTrack[j][i*fNstrips_ladder+k])-Dato; 

	DatoDec1=(UShort_t)(dato1*2);
	dato2=dato1*2-DatoDec1;
	
	DatoDec2=(UShort_t)(dato2*2);
	dato3=dato2*2-DatoDec2;
	
	DatoDec3=(UShort_t)(dato3*2);
	dato4=dato3*2-DatoDec3;
	
	DatoDec4=(UShort_t)(dato4*2);
	
	DatoDec=DatoDec1*0x0008+DatoDec2*0x0004+DatoDec3*0x0002+DatoDec4*0x0001;
      
	tempdat.push_back( ((Dato << 4) | (DatoDec & 0x000F)) );
	CkSum^=tempdat.back();
      };
      
      for (Int_t k=0; k<64;k++) { //0...63
        UShort_t DatoBad=0x0000;
        for (Int_t nb=0; nb<16;nb++) {
	  if( fSigmaTrack[j][i*fNstrips_ladder+k*16+nb]<0. ) DatoBad=( DatoBad | (0x8000 >> nb) ); 	 
	};

        tempdat.push_back(DatoBad);
	CkSum^=tempdat.back();
      };
      // end ladder

      // write in buffer the end ladder word

      switch(i){
      case 0: 
	tempdat.push_back(0x1807);
	break;
      case 1:
	tempdat.push_back(0x1808);
	break;	
      case 2:
	tempdat.push_back(0x1808);
	break;
      default:
	break;
      }
      CkSum^=tempdat.back();

      // write in buffer the TRAILER
      ReLength=(UShort_t)((fNstrips_ladder*2+64+1)*2+3);
      OveCheckCode=0x0000;
      

      tempdat.push_back(0x0000);
      tempdat.push_back((ReLength >> 8));
      tempdat.push_back(( (ReLength << 8) | (OveCheckCode & 0x00FF) ));
      
    };  // end TRAILER        
    
    cout<<"LENGTH OF TEMPDAT"<<tempdat.size()<<endl;
    // write in buffer the DSP header
    
    fData.push_back((0xE800 | ( ((j+1) << 3) | 0x0005) ));
    fData.push_back(0x01A9);
    fData.push_back(0x8740);
    EVENT_CAL=0;
    fData.push_back((0x1A00 | ( (0x03FF & EVENT_CAL)>> 1) ));
    PED_L1=0;
    fData.push_back(( ((EVENT_CAL << 15) | 0x5002 ) | ((0x03FF & PED_L1) << 2) ));
    fData.push_back(0x8014);
    fData.push_back(0x00A0);
    fData.push_back(0x0500);
    fData.push_back(0x2801);
    fData.push_back(0x400A);
    fData.push_back(0x0050);
    CkSum=(CkSum >> 8)^(CkSum&0x00FF);
    fData.push_back((0x0280 | (CkSum >> 3)));
    fData.push_back((0x1FFF | (CkSum << 13) ));
  
    ReLength=(UShort_t)((13*2)+3);
    OveCheckCode=0x0000;
    fData.push_back(0x0000);  
    fData.push_back((ReLength >> 8));
    fData.push_back(( (ReLength << 8) | (OveCheckCode & 0x00FF) ));
    cout<<"fDATA length before adding:"<<fData.size()<<endl;
    // Now we will copy to fData vector data from tempdat:
    p = tempdat.begin();
    while( p != tempdat.end() ){
      fData.push_back(*p);
      p++;
    }
    cout<<"fDATA length after adding:"<<fData.size()<<endl;
  };
}

  void PamVMCTrkDig::WriteCalib(){
    cout<<"Writing Tracker Calibrations..."<<endl;
    fraw->WritePSCU(&fDataPSCU);
    fraw->CopyUShortToBuff(&fData);
    if(fPadding) fraw->WritePadding(&fDataPadding);
}

void PamVMCTrkDig::LoadMipCor(){

std:: cout << "Entering LoadMipCor" << endl;

  Float_t xfactor=1./151.6*1.04;
  Float_t yfactor=1./152.1;

  fMipCor[0][0]=140.02*yfactor;
  fMipCor[0][1]=140.99*xfactor;
  fMipCor[0][2]=134.48*yfactor;
  fMipCor[0][3]=144.41*xfactor;
  fMipCor[0][4]=140.74*yfactor;
  fMipCor[0][5]=142.28*xfactor;
  fMipCor[0][6]=134.53*yfactor;
  fMipCor[0][7]=140.63*xfactor;
  fMipCor[0][8]=135.55*yfactor;
  fMipCor[0][9]=138.00*xfactor;
  fMipCor[0][10]=154.95*yfactor;
  fMipCor[0][11]=158.44*xfactor;
  
  
  fMipCor[1][0]=136.07*yfactor;
  fMipCor[1][1]=135.59*xfactor;
  fMipCor[1][2]=142.69*yfactor;
  fMipCor[1][3]=138.19*xfactor;
  fMipCor[1][4]=137.35*yfactor;
  fMipCor[1][5]=140.23*xfactor;
  fMipCor[1][6]=153.15*yfactor;
  fMipCor[1][7]=151.42*xfactor;
  fMipCor[1][8]=129.76*yfactor;
  fMipCor[1][9]=140.63*xfactor;
  fMipCor[1][10]=157.87*yfactor;
  fMipCor[1][11]=153.64*xfactor;

  fMipCor[2][0]=134.98*yfactor;
  fMipCor[2][1]=143.95*xfactor;
  fMipCor[2][2]=140.23*yfactor;
  fMipCor[2][3]=138.88*xfactor;
  fMipCor[2][4]=137.95*yfactor;
  fMipCor[2][5]=134.87*xfactor;
  fMipCor[2][6]=157.56*yfactor;
  fMipCor[2][7]=157.31*xfactor;
  fMipCor[2][8]=141.37*yfactor;
  fMipCor[2][9]=143.39*xfactor;
  fMipCor[2][10]=156.15*yfactor;
  fMipCor[2][11]=158.79*xfactor;

}


void PamVMCTrkDig::Digitize(){

  fData.clear();
 

  Int_t Iview;
  Int_t Nstrip;

  for (Int_t j=0; j<fNviews;j++) {

    for (Int_t i=0; i<fNladder;i++) {

      Float_t commonN1=gRandom->Gaus(0.,fSigmaCommon);
      Float_t commonN2=gRandom->Gaus(0.,fSigmaCommon);
      for (Int_t k=0; k<fNstrips_ladder;k++) {
      Nstrip=i*fNstrips_ladder+k;
      Float_t Sigma=Abs(fSigmaTrack[j][Nstrip]);
      AdcTrack[j][Nstrip]=gRandom->Gaus(fPedeTrack[j][Nstrip],Sigma );
      if(k<4*128) {AdcTrack[j][Nstrip] += commonN1;}  // full correlation of 4 VA1 Com. Noise
      else {AdcTrack[j][Nstrip] += commonN2;}   // full correlation of 4 VA1 Com. Noise
      if(AdcTrack[j][Nstrip] < 0. ) AdcTrack[j][Nstrip]=0.;
      if(AdcTrack[j][Nstrip] > 4095.) AdcTrack[j][Nstrip]=4095.;
      };  
    };
  }; 

  if (fhits){
    Float_t ADCfull;
    Int_t iladd=0;
    for (Int_t ix=0; ix<fhits->GetNXHit();ix++) {
      Iview=((fhits->GetXHit(ix))->fnpstrip)*2-1;
      Nstrip=(Int_t)((fhits->GetXHit(ix))->fistrip)-1; 
      if(Nstrip<fNstrips_ladder) iladd=0;
      if((Nstrip>=fNstrips_ladder)&&(Nstrip<2*fNstrips_ladder)) iladd=1;
      if((Nstrip>=2*fNstrips_ladder)&&(Nstrip<3*fNstrips_ladder)) iladd=2;
      ADCfull=AdcTrack[Iview][Nstrip] += ((fhits->GetXHit(ix))->fqstrip)*fMipCor[iladd][Iview];
      AdcTrack[Iview][Nstrip] *= SaturationTrackx(ADCfull);
    };
    
    
    for (Int_t iy=0; iy<fhits->GetNYHit();iy++) {
      Iview=((fhits->GetYHit(iy))->fnpstrip)*2-2;
      Nstrip=(Int_t)((fhits->GetYHit(iy))->fistrip)-1; 
      if(Nstrip<fNstrips_ladder) iladd=0;
      if((Nstrip>=fNstrips_ladder)&&(Nstrip<2*fNstrips_ladder)) iladd=1;
      if((Nstrip>=2*fNstrips_ladder)&&(Nstrip<3*fNstrips_ladder)) iladd=2;
      ADCfull=AdcTrack[Iview][Nstrip] -= ((fhits->GetYHit(iy))->fqstrip)*fMipCor[iladd][Iview];
      AdcTrack[Iview][Nstrip] *= SaturationTracky(ADCfull);
    };

    CompressTrackData();
  };
}


void PamVMCTrkDig::CompressTrackData(){

// copy of the corresponding compression fortran routine + new digitization

  Int_t oldval=0;
  Int_t newval=0;
  Int_t trasmesso=0;
  Int_t ntrastot=0;
  Float_t real, inte;
  Int_t cercacluster=0;
  Int_t kt=0;
  static const int DSPbufferSize = 4000; // 13 bit buffer to be rearranged in 16 bit Track buffer
  UShort_t DataDSP[DSPbufferSize];  // 13 bit  buffer to be rearranged in 16 bit Track buffer
  UShort_t DSPlength;  // 13 bit buffer to be rearranged in 16 bit Track buffer
  UShort_t CheckSum, Nword, Dato, DATA, ReLength, OveCheckCode; 
  Int_t k, diff, clval, clvalp, klp, kl1, kl2, Bit16free, Bit13ToWrite;


  for (Int_t iv=0; iv<fNviews;iv++) {
    memset(DataDSP,0,sizeof(UShort_t)*DSPbufferSize);
    DSPlength=16; // skip the header, to be written later
    CheckSum=0;

    for (Int_t ladder=0; ladder<fNladder;ladder++) {
      k=0;
      while (k<fNstrips_ladder) {
	// compress write in buffer the current LADDER
	if ( k == 0)  {
	  real=modff(AdcTrack[iv][ladder*fNstrips_ladder+k],&inte);
	  if (real > 0.5) inte=inte+1;
	  newval=(Int_t)inte -(Int_t)fPedeTrack[iv][ladder*fNstrips_ladder+k];
	  // first strip of ladder is transmitted
	  DataDSP[DSPlength]=( ((UShort_t)inte) & 0x0FFF);
	  DSPlength++;  
	  ntrastot++;
	  trasmesso=1;
	  oldval=newval;
	  kt=k;
	  k++;
	  continue; 
	};
	real=modff(AdcTrack[iv][ladder*fNstrips_ladder+k],&inte);
	if (real > 0.5) inte=inte+1;
	newval=(Int_t)inte -(Int_t)(fPedeTrack[iv][ladder*fNstrips_ladder+k]);
	
	cercacluster=1; //PAR? 
	
	if (cercacluster==1) {
	  
	  diff=0;	  
	  switch ((iv+1)%2) {
	  case 0: diff=newval-oldval;
	    break;
	  case 1: diff=oldval-newval;	  
	    break;
	  };
	  
	  
	  if (diff>fCutclu*(Int_t)fabs(fSigmaTrack[iv][ladder*fNstrips_ladder+k]) ) {
	    clval=newval;
	    klp=k;	// go on to search for maximum 
	    klp++;
	    
	    while(klp<fNstrips_ladder) {
	      real=modff(AdcTrack[iv][ladder*fNstrips_ladder+klp],&inte);
	      if (real > 0.5) inte=inte+1;
	      clvalp=(Int_t)inte -(Int_t)fPedeTrack[iv][ladder*fNstrips_ladder+klp];
	      if((iv+1)%2==0) {
		
		if(clvalp>clval) {
		  clval=clvalp;
		  k=klp;}
		else break; // max of cluster found
		
	      } else {  
		
		if(clvalp<clval) {
		  clval=clvalp;
		  k=klp;}
		else break; // max of cluster found       
	      };
	      klp++;
	    };
	    
            kl1=k-fNclst; // max of cluster (or end of ladder ?)
	    trasmesso=0;
	    if(kl1<0)  kl1=0;
	    if(kt>=kl1) kl1=kt+1;
            if( (kt+1)==kl1 ) trasmesso=1;
	    
	    kl2=k+fNclst;
	    if(kl2>=fNstrips_ladder) kl2=fNstrips_ladder-1;
	    
	    for(Int_t klt=kl1 ; klt<=kl2 ; klt++) {
              if(trasmesso==0) { 
		
	        DataDSP[DSPlength]=( ((UShort_t)klt) | 0x1000);
	        DSPlength++;  
                ntrastot++;
		
	        real=modff(AdcTrack[iv][ladder*fNstrips_ladder+klt],&inte);
                if (real > 0.5) inte=inte+1;
	        DataDSP[DSPlength]=( ((UShort_t)inte) & 0x0FFF);
	        DSPlength++; 
	        ntrastot++;
		
	      }
	      else { 
		
 	        real=modff(AdcTrack[iv][ladder*fNstrips_ladder+klt],&inte);
                if (real > 0.5) inte=inte+1;
	        DataDSP[DSPlength]=( ((UShort_t)inte) & 0x0FFF);
	        DSPlength++;               
		ntrastot++;
	      };
	      trasmesso=1;	      	   	    
	    };  // end trasmission
	    kt=kl2;
	    k=kl2;
	    real=modff(AdcTrack[iv][ladder*fNstrips_ladder+kt],&inte);
            if (real > 0.5) inte=inte+1;
 	    oldval=(Int_t)inte -(Int_t)fPedeTrack[iv][ladder*fNstrips_ladder+kt];
            k++;
	    continue;
	 }
       }//END cercacluster

       // start ZOP check for strips no

       if(abs(newval-oldval)>=fCutzop*(Int_t)fabs(fSigmaTrack[iv][ladder*fNstrips_ladder+k]) ) {	 if(trasmesso==0) {
	 	 
	 DataDSP[DSPlength]=( ((UShort_t)k) | 0x1000);
	 DSPlength++;  
         ntrastot++;    
	 
	 real=modff(AdcTrack[iv][ladder*fNstrips_ladder+k],&inte);
         if (real > 0.5) inte=inte+1;
	 DataDSP[DSPlength]=( ((UShort_t)inte) & 0x0FFF);
	 DSPlength++; 
	 ntrastot++;

       } else {
	 real=modff(AdcTrack[iv][ladder*fNstrips_ladder+k],&inte);
         if (real > 0.5) inte=inte+1;
	 DataDSP[DSPlength]=(  ((UShort_t)inte) & 0x0FFF);
	 DSPlength++; 
	 ntrastot++;      
       };
       trasmesso=1;
       oldval=newval;
       kt=k;
       
       }  
       else trasmesso=0;
       // end zop
       k++;      
     };

     DataDSP[DSPlength]=( ((UShort_t)(ladder+1))  | 0x1800);
     DSPlength++; 
     ntrastot++; 
     trasmesso=0;
     
   };  //end cycle inside dsp 

   // here put DSP header
   DataDSP[0]=(0x1CA0 | ((UShort_t)(iv+1)) );
   Nword=(DSPlength*13)/16; 
   if( ((DSPlength*13)%16)!=0) Nword++;
   DataDSP[1]=(0x1400 | ( Nword >> 10));
   DataDSP[2]=(0x1400 | ( Nword & 0x03FF) );
   DataDSP[3]=(0x1400 | (( (UShort_t)(fraw->GetCounter() >> 10) ) & 0x03FF) ); 
   DataDSP[4]=(0x1400 | (( (UShort_t)(fraw->GetCounter()) )  & 0x03FF) );
   DataDSP[5]=(0x1400 | ( (UShort_t)(fNclst << 7) ) | ( (UShort_t)(fCutzop << 4) )
	       | ( (UShort_t)fCutzop  ) ); 
   DataDSP[6]=0x1400;
   DataDSP[7]=0x1400;
   DataDSP[8]=0x1400;
   DataDSP[9]=0x1400;
   DataDSP[10]=0x1400;
   DataDSP[11]=0x1400;
   DataDSP[12]=0x1400;
   DataDSP[13]=0x1400;
   DataDSP[14]=(0x1400 | (CheckSum & 0x00FF) ); 
   DataDSP[15]=0x1C00;
   // end DSP header    

   // write 13 bit DataDSP bufer inside 16 bit fDataTrack buffer

   Bit16free=16;
   DATA = 0;
   for (Int_t NDSP=0; NDSP<DSPlength;NDSP++) {
     Bit13ToWrite=13;
     while(Bit13ToWrite>0) { 
       if(Bit13ToWrite<=Bit16free) {
	 Dato=((DataDSP[NDSP]&(0xFFFF >> (16-Bit13ToWrite)))<<(Bit16free-Bit13ToWrite));
	 DATA = DATA | Dato ;
	 Bit16free=Bit16free-Bit13ToWrite;
	 Bit13ToWrite=0;
	 if(Bit16free==0) {
	   if(NDSP>15) CheckSum=CheckSum^DATA;
	   fData.push_back(DATA);
	   DATA = 0;
	   Bit16free=16;
	 };          
       }
       else if(Bit13ToWrite>Bit16free) {
	 Dato=( (DataDSP[NDSP]&(0xFFFF >> (16-Bit13ToWrite) ) ) >> (Bit13ToWrite-Bit16free) );
	 DATA = DATA | Dato ;
	 fData.push_back(DATA);
	 if(NDSP>15) CheckSum=CheckSum^DATA;
	 DATA = 0;
	 Bit13ToWrite=Bit13ToWrite-Bit16free;
	 Bit16free=16;	 
       };
       
     }; // end cycle while(Bit13ToWrite>0)
     
   }; // end cycle DataDSP

   if(Bit16free!=16) { 
     fData.push_back(DATA); 
     DATA= 0; 
     CheckSum=CheckSum^DATA; 
   };
   CheckSum=(CheckSum >> 8)^(CheckSum&0x00FF); 
   fData.at(fData.size()-Nword+11)=(0x0280 | (CheckSum >> 3));
   fData.at(fData.size()-Nword+12)=(0x1C00 | (CheckSum << 13) );
   
   // end write 13 bit DataDSP bufer inside 16 bit fDataTrack buffer

   //write trailer on buffer
   ReLength=(UShort_t)((Nword+13)*2+3);
   OveCheckCode=0x0000;
   fData.push_back(0x0000);
   fData.push_back((ReLength >> 8));
   fData.push_back(( (ReLength << 8) | (OveCheckCode & 0x00FF) ));
   
   // end trailer 


 }//END VIEWS CYCLE


}


Float_t PamVMCTrkDig::SaturationTrackx(Float_t ADC) {
  Float_t SatFact=1.;
  if(ADC<1.) { SatFact=1./ADC; }; 
  if(ADC>3000.) { SatFact=3000./ADC; }; 
  return SatFact;
};


Float_t PamVMCTrkDig::SaturationTracky(Float_t ADC) {
  Float_t SatFact=1.;
  if(ADC<70.) { SatFact=70./ADC; }; 
  if(ADC>4095.) { SatFact=4095./ADC; }; 
  return SatFact;
};