source: trunk/MagicSoft/Simulation/Detector/StarResponse/MStarLight.cxx@ 5099

#include "MStarLight.hxx" 

using namespace std;

MStarLight::MStarLight() {
  //
  //   default constructor
  //

  fBrightness  = 0. ; 
  fTimeRange   = TIMERANGE  ; 

  fBinsTrig    = TRIGBINS   ; 
  fTrigShape   = 0. ; 
  fAmplTrig    = 0. ; 
  fFwhmTrig    = 0. ; 

  fBinsFadc    = FADCBINS   ; 
  fFadcShape   = 0. ; 
  fIntegFadc    = 0. ; 
  fFwhmFadc    = 0. ; 

  
  for (Int_t i= 0; i< TRIGBINS ; i++) 
    fTrig[i] = 0. ; 

  for (Int_t i= 0; i<FADCBINS; i++) 
    fFadc[i] = 0. ; 
  
  for (Int_t i= 0; i<40; i++) 
    fTrigResp[i] = 0. ; 

  for (Int_t i= 0; i<45; i++) 
    fFadcResp[i] = 0. ; 
}

void MStarLight::Reset() {

  
  fBrightness  = 0. ; 
  fTimeRange   = TIMERANGE  ; 

  fBinsTrig    = TRIGBINS   ; 
  fTrigShape   = 0. ; 
  fAmplTrig    = 0. ; 
  fFwhmTrig    = 0. ; 

  fBinsFadc    = FADCBINS   ; 
  fFadcShape   = 0. ; 
  fIntegFadc    = 0. ; 
  fFwhmFadc    = 0. ; 

  
  for (Int_t i= 0; i< TRIGBINS ; i++) 
    fTrig[i] = 0. ; 

  for (Int_t i= 0; i<FADCBINS; i++) 
    fFadc[i] = 0. ; 
  
  for (Int_t i= 0; i<40; i++) 
    fTrigResp[i] = 0. ; 

  for (Int_t i= 0; i<45; i++) 
    fFadcResp[i] = 0. ; 
}



Float_t MStarLight::GetBrightness () 
{
  return fBrightness  ; 
}

void MStarLight::SetBrightness (Float_t in ) 
{
  fBrightness = in ; 
}

Float_t MStarLight::GetAmplTrig () 
{
  return fAmplTrig  ; 
}

void MStarLight::SetAmplTrig (Float_t in ) 
{
  fAmplTrig = in ; 
}

Float_t  MStarLight::GetFwhmTrig () 
{
  return fFwhmTrig  ; 
}

void MStarLight::SetFwhmTrig (Float_t in ) 
{
  fFwhmTrig = in ; 
}


Float_t MStarLight::GetIntegFadc () 
{
  return fIntegFadc  ; 
}

void MStarLight::SetIntegFadc (Float_t in ) 
{
  fIntegFadc = in ; 
}

Float_t  MStarLight::GetFwhmFadc () 
{
  return fFwhmFadc  ; 
}

void MStarLight::SetFwhmFadc (Float_t in ) 
{
  fFwhmFadc = in ; 
}


void MStarLight::SetTrigResponse( Float_t *in ) {
  
  for (Int_t i= 0; i<40; i++) 
    fTrigResp[i] = in[i] ; 
  
}

void MStarLight::SetFadcResponse( Float_t *in ) {

  for (Int_t i= 0; i<45; i++)
    fFadcResp[i] = in[i] ; 
}

void MStarLight::FillResponse( Float_t ampl, Float_t time ) {

  //     fill the trigger response
  
  Int_t startbin = (Int_t) (time * 4)  ; 
  
  Int_t icount = 0 ;
  Int_t idata      ; 
  
  for ( Int_t i=startbin ; i<startbin+40 ; i++) {
    //
    if ( i < TRIGBINS ) 
      idata = i ; 
    else if ( i >= TRIGBINS ) 
      idata = i - TRIGBINS ; 

    fTrig[idata] = fTrig[idata] + ampl * fTrigResp[icount] ; 
    
    icount++ ; 
  }

  //
  //   fill the FADC content
  //

  startbin = (Int_t) ( time / 0.6666666666 );
  
  Int_t ichanfadc = 0 ; 
  
  //
  //   putting the response slices in the right sig slices.
  //   Be carefull, because both slices have different widths. 
  //

  //
  // Changed, Jan 2004, A. Moralejo: now the FADC does not integrate all
  // the signal within each FADC slice, but measures just the signal height
  // at one point, like the real FADC does. By default, each FADC slice 
  // contains 5 bins of the response histogram(fFadcResp[45]). Warning:
  // do not change this unless you do the corresponding modifications also
  // in MFadc.cxx, or the signal and the noise photoelectrons will be treated 
  // differently!!
  //

  for ( Int_t i = 0 ; i < 45; i += 5 ) {
    ichanfadc = (Int_t) ((startbin+i)/5) ; 
    
    if ( (ichanfadc) < FADCBINS ) 
      idata = ichanfadc ; 
    else if ( ichanfadc >=FADCBINS ) 
      idata = ichanfadc-FADCBINS ; 
    
    fFadc[idata] += (ampl * fFadcResp[i] )  ; 
    
  }

}

void MStarLight::ElecNoise ( Float_t noiseTrig, Float_t noiseFadc ) {
  //
  //     putting some noise to the baseline
  // 

  TRandom2 wuerfel( (UInt_t) (noiseTrig*100) ) ; 
 
  for (Int_t i=0; i< TRIGBINS ; i++ ) {
    fTrig[i] += wuerfel.Gaus(0., noiseTrig ) ; 
  }

  for (Int_t i=0; i< FADCBINS ; i++ ) {
    fFadc[i] += wuerfel.Gaus(0., noiseFadc ) ; 
  }
  
  
}

Float_t MStarLight::GetTrig( Int_t i){
  //------------------------------------------------------------------
  //
  //  It gets the value of the simualted trigger in the i bin
  //
  return fTrig[i];
}

Float_t MStarLight::GetFadc( Int_t i){
  //------------------------------------------------------------------
  //
  //  It gets the value of the simualted FADC signal in the i bin
  //
  return fFadc[i];
}

void MStarLight::StoreHisto( char *filename) {
  
  Float_t baseline = 0. ; 


  //    first the histogramms for trigger
  //
  //    the calculated trigger signal before baseline
  // 

  TH1F  trigresp ("trigresp", "Trigger Response", TRIGBINS, 0., TIMERANGE) ; 

  for (Int_t i=0; i< TRIGBINS ; i++ ) {
    
    trigresp.SetBinContent(i, fTrig[i]) ; 

    baseline = baseline + fTrig[i]  ; 
  }

  baseline = baseline / TRIGBINS ; 
  
  TH1F  trigbase ("trigbase", "Response after Baseline shift", 
                  TRIGBINS, 0., TIMERANGE) ; 
  for (Int_t i=0; i< TRIGBINS ; i++ ) {
    
    trigbase.SetBinContent(i, fTrig[i]-baseline ) ; 
    
  }
  
  TH1F  trigdist ("trigdist", "Noise on the baseline", 
                  1000, -25., 25. ) ; 

  for (Int_t i=0; i< TRIGBINS ; i++ ) {
    
    trigdist.Fill( (Float_t) trigbase.GetBinContent(i) ) ; 
  }


  // 
  //     second the histograms for the fadc
  //
  
  
  TH1F  fadcresp ("fadcresp", "Fadc Response", FADCBINS, 0., TIMERANGE) ; 

  for (Int_t i=0; i< FADCBINS ; i++ ) {
    
    fadcresp.SetBinContent(i, fFadc[i]) ; 

    baseline = baseline + fFadc[i]  ; 
  }

  baseline = baseline / FADCBINS ; 
  
  TH1F  fadcbase ("fadcbase", "Fadc after Baseline shift", 
                  FADCBINS, 0., TIMERANGE) ; 
  for (Int_t i=0; i< FADCBINS ; i++ ) {
    
    fadcbase.SetBinContent(i, fFadc[i]-baseline ) ; 
    
  }
  
  TH1F  fadcdist ("fadcdist", "Noise on fadc", 
                  1000, -100., 100. ) ; 

  for (Int_t i=0; i< FADCBINS ; i++ ) {
    
    fadcdist.Fill( (Float_t) fadcbase.GetBinContent(i) ) ; 
  }



  
  TFile outfile( filename, "UPDATE" ) ; 
  
  trigresp.Write() ; 
  trigbase.Write() ; 
  trigdist.Write() ; 

  fadcresp.Write() ; 
  fadcbase.Write() ; 
  fadcdist.Write() ; 
  
  
  outfile.Close() ; 
}

void MStarLight::WriteBinary( char *filename) {
  //
  //  write the things to the binary file

  FILE *datei ; 

  datei = fopen ( filename, "w" ) ; 
  
  if ( ! datei ) {
    cout << " ERROR: Can't open the file " << filename 
         << endl ; 
    exit (230) ; 
  }

  Float_t version = VERSIONSR;

  //   write them out 

  fwrite ( &version, sizeof(Float_t), 1, datei ) ; 
  fwrite ( &fBrightness, sizeof(Float_t), 1, datei ) ; 
  fwrite ( &fTimeRange , sizeof(Float_t), 1, datei ) ; 

  fwrite ( &fBinsTrig  , sizeof(Float_t), 1, datei ) ; 
  fwrite ( &fTrigShape , sizeof(Float_t), 1, datei ) ; 
  fwrite ( &fAmplTrig  , sizeof(Float_t), 1, datei ) ; 
  fwrite ( &fFwhmTrig  , sizeof(Float_t), 1, datei ) ; 

  fwrite ( &fBinsFadc  , sizeof(Float_t), 1, datei ) ; 
  fwrite ( &fFadcShape , sizeof(Float_t), 1, datei ) ; 
  fwrite ( &fIntegFadc , sizeof(Float_t), 1, datei ) ; 
  fwrite ( &fFwhmFadc  , sizeof(Float_t), 1, datei ) ; 

  fwrite ( &fTrig      , sizeof(Float_t), TRIGBINS, datei ) ; 

  //    We want to store the FADC signal taking into account the AC
  //    coupling
  //
  //    We calculate and substract the baseline
  // 

  Float_t baseline = 0. ; 

  for (Int_t i=0; i< FADCBINS ; i++ ) {
    
    baseline = baseline + fFadc[i]  ; 
  }

  baseline = baseline / FADCBINS ; 
  
  for (Int_t i=0; i< FADCBINS ; i++ ) {
    
      fFadc[i]=fFadc[i]-baseline;
    
  }
  
  fwrite ( &fFadc      , sizeof(Float_t), FADCBINS, datei ) ; 


  fclose ( datei ) ; 
  
}

void MStarLight::ReadBinary( char *filename) {
  
  //
  //  read the things from the binary file

  FILE *datei ; 

  datei = fopen ( filename, "r" ) ; 
  
  if ( ! datei ) {
    cout << " ERROR: Can't open the file " << filename 
         << endl ; 
    cout<<  "        The database for the NSB may be too small. "
        <<  endl;
    cout<<  "        See the How to Use of the Camera simulation for more information"
        <<  endl;
    exit (230) ; 
  }

  Float_t read_version;
  Float_t current_version=VERSIONSR;

  // Check taht we read the right version of the Database
  fread ( &read_version, sizeof(Float_t), 1, datei ) ; 

  if(Int_t(read_version)!=Int_t(current_version)){
    cout<<" ERROR: You are trying to read database VERSION "<<
      read_version<<endl;
    cout<<"        You must generate a database for the current VERSION "<<
      current_version<<endl;
    cout<<"        See the NSB database section in the Camera How to Use note."
        <<endl;
    exit (230);
  }

  fread ( &fBrightness, sizeof(Float_t), 1, datei ) ; 
  fread ( &fTimeRange , sizeof(Float_t), 1, datei ) ; 


  fread ( &fBinsTrig  , sizeof(Float_t), 1, datei ) ; 
  fread ( &fTrigShape , sizeof(Float_t), 1, datei ) ; 
  fread ( &fAmplTrig  , sizeof(Float_t), 1, datei ) ; 
  fread ( &fFwhmTrig  , sizeof(Float_t), 1, datei ) ; 

  fread ( &fBinsFadc  , sizeof(Float_t), 1, datei ) ; 
  fread ( &fFadcShape , sizeof(Float_t), 1, datei ) ; 
  fread ( &fIntegFadc , sizeof(Float_t), 1, datei ) ; 
  fread ( &fFwhmFadc  , sizeof(Float_t), 1, datei ) ; 

  fread ( &fTrig      , sizeof(Float_t), TRIGBINS, datei ) ; 

  fread ( &fFadc      , sizeof(Float_t), FADCBINS, datei ) ; 

  

  fclose ( datei ) ; 
}