]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - START/AliSTARTDigitizer.cxx
- The part of JETAN dealing with ESD data has been separated from the one using MC...
[u/mrichter/AliRoot.git] / START / AliSTARTDigitizer.cxx
index 25caf7e20f355e6aa9fcd97c12d4e6b418b05ce3..faf408424ca8470a47a511c3c789bdcb2a7a2b3c 100644 (file)
@@ -22,6 +22,7 @@
 #include <TArrayI.h>
 #include <TError.h>
 #include <TH1F.h>
+#include <TGraph.h>
 
 #include "AliLog.h"
 #include "AliSTARTDigitizer.h"
 #include <stdlib.h>
 #include <Riostream.h>
 #include <Riostream.h>
+#include "AliSTARTParameters.h"
+#include "AliCDBLocal.h"
+#include "AliCDBStorage.h"
+#include "AliCDBManager.h"
+#include "AliCDBEntry.h"
 
 ClassImp(AliSTARTDigitizer)
 
@@ -52,14 +58,11 @@ AliSTARTDigitizer::AliSTARTDigitizer(AliRunDigitizer* manager)
    fSTART(0),
    fHits(0),
    fdigits(0),
-   ftimeTDC(0),
+   ftimeCFD(0),
+   ftimeLED(0),
    fADC(0),
-   ftimeTDCAmp(0),
-   fADCAmp(0),
-   fSumMult(0),
-   fEff(0)
+   fADC0(0)
 {
-  //   cout<<"AliSTARTDigitizer::AliSTARTDigitizer"<<endl;
 // ctor which should be used
 
   AliDebug(1,"processed");
@@ -68,18 +71,12 @@ AliSTARTDigitizer::AliSTARTDigitizer(AliRunDigitizer* manager)
   fHits = 0;
   fdigits = 0;
 
-  ftimeTDC = new TArrayI(24); 
+  ftimeCFD = new TArrayI(24); 
   fADC = new TArrayI(24); 
-  ftimeTDCAmp = new TArrayI(24); 
-  fADCAmp = new TArrayI(24); 
-  fSumMult = new TArrayI(6); 
+  ftimeLED = new TArrayI(24); 
+  fADC0 = new TArrayI(24); 
   
 
-  TFile* file = TFile::Open("$ALICE_ROOT/START/PMTefficiency.root");
-  fEff = (TH1F*) file->Get("hEff")->Clone();
-  fEff->SetDirectory(NULL);
-  file->Close();
-  delete file;
 }
 
 //------------------------------------------------------------------------
@@ -89,12 +86,10 @@ AliSTARTDigitizer::~AliSTARTDigitizer()
 
   AliDebug(1,"START");
 
-  delete ftimeTDC;
+  delete ftimeCFD;
   delete fADC;
-  delete fEff;
-  delete ftimeTDCAmp;
-  delete  fADCAmp;
-  delete fSumMult;
+  delete ftimeLED;
+  delete  fADC0;
 }
 
  //------------------------------------------------------------------------
@@ -124,49 +119,78 @@ void AliSTARTDigitizer::Exec(Option_t* /*option*/)
   //output loader 
   AliRunLoader *outRL = AliRunLoader::GetRunLoader(fManager->GetOutputFolderName());
   AliLoader * pOutStartLoader = outRL->GetLoader("STARTLoader");
-  AliSTART *fSTART  = (AliSTART*)outRL ->GetAliRun()->GetDetector("START");
 
   AliDebug(1,"start...");
-  cout<<" AliSTARTDigitizer::Exec "<<endl;
   //input loader
   //
   // From hits to digits
   //
 Int_t hit, nhits;
-  Int_t countE[24], sumMult[3];
-  Int_t volume,pmt,tr,qt,qtAmp;
-  Int_t  bestRightTDC,bestLeftTDC;
-  Float_t time[24],besttime[24], timeGaus[24] ;
-  Float_t channelWidth=25.; //ps
+ Int_t hit, nhits;
+  Int_t countE[24];
+  Int_t volume, pmt, trCFD, trLED; 
+  Float_t sl, qt;
+  Int_t  bestRightTDC, bestLeftTDC, qtCh;
+  Float_t time[24], besttime[24], timeGaus[24] ;
     //Q->T-> coefficients !!!! should be asked!!!
-  Float_t ph2mV = 150./500.;
+  Float_t gain[24],timeDelayCFD[24], timeDelayLED[24];
   Int_t threshold =50; //photoelectrons
-  Int_t mV2channel=200000/(25*25);  //5V -> 200ns
+  Float_t zdetA, zdetC;
+   Int_t sumMultCoeff = 100;
+  TObjArray slewingLED;
+  TObjArray slewingRec;
+  AliSTARTParameters* param = AliSTARTParameters::Instance();
+  param->Init();
 
+  Int_t ph2Mip = param->GetPh2Mip();     
+  Int_t channelWidth = param->GetChannelWidth() ;  
+  Float_t delayVertex = param->GetTimeDelayTVD();
+  for (Int_t i=0; i<24; i++){
+    timeDelayCFD[i] = param->GetTimeDelayCFD(i);
+    timeDelayLED[i] = param->GetTimeDelayLED(i);
+    gain[i] = param->GetGain(i);
+    TGraph* gr = param ->GetSlew(i);
+    slewingLED.AddAtAndExpand(gr,i);
+
+    TGraph* gr1 = param ->GetSlewRec(i);
+    slewingRec.AddAtAndExpand(gr1,i);
+
+    TGraph* grEff = param ->GetPMTeff(i);
+    fEffPMT.AddAtAndExpand(grEff,i);
+  }
+  zdetC = param->GetZposition(0);
+  zdetA = param->GetZposition(1);
   
   AliSTARThit  *startHit;
   TBranch *brHits=0;
-
+  
   Int_t nFiles=fManager->GetNinputs();
   for (Int_t inputFile=0; inputFile<nFiles;  inputFile++) {
     if (inputFile < nFiles-1) {
-      AliWarning(Form("ignoring input stream %d", inputFile));
-      continue;
-      
+       AliWarning(Form("ignoring input stream %d", inputFile));
+       continue;
+       
     }
     
     Float_t besttimeright=99999.;
     Float_t besttimeleft=99999.;
-    Int_t timeDiff, meanTime;
-    
+    Int_t pmtBestRight=9999;
+    Int_t pmtBestLeft=9999;
+    Int_t timeDiff=999, meanTime=0;
+    Int_t sumMult =0;   fSumMult=0;
+    bestRightTDC = 99999;  bestLeftTDC = 99999;
+    ftimeCFD -> Reset();
+    fADC -> Reset();
+    fADC0 -> Reset();
+    ftimeLED ->Reset();
     for (Int_t i0=0; i0<24; i0++)
       {
-       time[i0]=99999;  besttime[i0]=99999;    countE[i0]=0;
+       time[i0]=besttime[i0]=timeGaus[i0]=99999; countE[i0]=0;
       }
-    for ( Int_t imu=0; imu<3; imu++) sumMult[imu]=0;
     AliRunLoader * inRL = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(inputFile));
     AliLoader * pInStartLoader = inRL->GetLoader("STARTLoader");
     if (!inRL->GetAliRun()) inRL->LoadgAlice();
+    AliSTART *fSTART  = (AliSTART*)inRL ->GetAliRun()->GetDetector("START");
 
        //read Hits 
     pInStartLoader->LoadHits("READ");//probably it is necessary to load them before
@@ -195,82 +219,102 @@ void AliSTARTDigitizer::Exec(Option_t* /*option*/)
          }
          pmt=startHit->Pmt();
          Int_t numpmt=pmt-1;
-                 Double_t e=startHit->Etot();
+         Double_t e=startHit->Etot();
          volume = startHit->Volume();
          
-         if(e>0 && RegisterPhotoE(e)) {
+         //      if(e>0 && RegisterPhotoE(numpmt,e)) {
+         if(e>0 ) {
            countE[numpmt]++;
            besttime[numpmt] = startHit->Time();
            if(besttime[numpmt]<time[numpmt])
              {
                time[numpmt]=besttime[numpmt];
              }
-         } 
+         } //photoelectron accept 
        } //hits loop
     } //track loop
     
     //spread time right&left by 25ps   && besttime
+    Float_t c = 0.0299792; // cm/ps
+    
+    Float_t koef=(zdetA-zdetC)/c; //correction position difference by cable
     for (Int_t ipmt=0; ipmt<12; ipmt++){
       if(countE[ipmt] > threshold) {
-       timeGaus[ipmt]=gRandom->Gaus(time[ipmt],25);
-       if(timeGaus[ipmt]<besttimeleft) besttimeleft=timeGaus[ipmt]; //timeleft
-       sumMult[0] +=  countE[ipmt];
-       sumMult[1] +=  countE[ipmt];
-      }
+       timeGaus[ipmt]=gRandom->Gaus(time[ipmt],25)+koef;
+       if(timeGaus[ipmt]<besttimeleft){
+         besttimeleft=timeGaus[ipmt]; //timeleft
+         pmtBestLeft=ipmt;}
+     }
     }
-    for ( Int_t ipmt=12; ipmt<24; ipmt++){
+     for ( Int_t ipmt=12; ipmt<24; ipmt++){
       if(countE[ipmt] > threshold) {
        timeGaus[ipmt]=gRandom->Gaus(time[ipmt],25); 
-       if(timeGaus[ipmt]<besttimeright)  besttimeright=timeGaus[ipmt]; //timeright
-       sumMult[0] +=  countE[ipmt];
-       sumMult[2] +=  countE[ipmt];
+       if(timeGaus[ipmt]<besttimeright) {
+         besttimeright=timeGaus[ipmt]; //timeright
+         pmtBestRight=ipmt;}
       }        
     }
-    //ADC features fill digits
-   //folding with experimental time distribution
-    Float_t c = 0.0299792; // cm/ps
-    Float_t koef=(350.-69.7)/c;
-    Float_t  besttimeleftR= besttimeleft;
-    besttimeleft=koef+besttimeleft;
-    bestLeftTDC=Int_t (besttimeleftR/channelWidth);
-    bestRightTDC=Int_t (besttimeright/channelWidth);
-    timeDiff=Int_t ((besttimeright-besttimeleftR)/channelWidth);
-    meanTime=Int_t (((besttimeright+besttimeleft)/2.)/channelWidth);
-    AliDebug(2,Form(" time in ns %f ", besttimeleft));
+   //folding with alignmentz position distribution  
+    if( besttimeleft > 10000. && besttimeleft <15000)
+      bestLeftTDC=Int_t ((besttimeleft+1000*timeDelayCFD[pmtBestLeft])
+                        /channelWidth);
+    if( besttimeright > 10000. && besttimeright <15000)
+      bestRightTDC=Int_t ((besttimeright+1000*timeDelayCFD[pmtBestRight])
+                       /channelWidth);
+
+    if (bestRightTDC < 99999 && bestLeftTDC < 99999)
+      {
+       timeDiff=Int_t (((besttimeleft-besttimeright)+1000*delayVertex)
+                       /channelWidth);
+       meanTime=Int_t (((besttimeright+1000*timeDelayCFD[pmtBestLeft]+
+                         besttimeleft+1000*timeDelayCFD[pmtBestLeft])/2.)
+                       /channelWidth);
+      }
+       AliDebug(10,Form(" time right& left %i %i  time diff && mean time in channels %i %i",bestRightTDC,bestLeftTDC, timeDiff, meanTime));
     for (Int_t i=0; i<24; i++)
       {
-       //  fill TDC
-       tr= Int_t (timeGaus[i]/channelWidth); 
-       if(timeGaus[i]>50000) tr=0;
-       
-       //fill ADC
-               Int_t al= countE[i]; 
-       // QTC procedure:
-       // phe -> mV 0.3; 1MIP ->500phe -> ln (amp (mV)) = 5;
-       // max 200ns, HIJING  mean 50000phe -> 15000mv -> ln = 15 (s zapasom)
-       // channel 25ps
-       if (al>threshold) {
-         qt=Int_t (TMath::Log(al*ph2mV) * mV2channel); 
-         qtAmp=Int_t (TMath::Log(al*10*ph2mV) * mV2channel);
-         cout<<i<<" "<<qt<<" "<<qtAmp<<endl;
-         fADC->AddAt(qt,i);
-         ftimeTDC->AddAt(tr,i);
-         fADCAmp->AddAt(qtAmp,i);
-         ftimeTDCAmp->AddAt(tr,i); //now is the same as non-amplified but will be change
+               Float_t  al = countE[i]; 
+       if (al>threshold && timeGaus[i]<50000 ) {
+         //fill ADC
+         // QTC procedure:
+         // phe -> mV 0.3; 1MIP ->500phe -> ln (amp (mV)) = 5;
+         // max 200ns, HIJING  mean 50000phe -> 15000mv -> ln = 15 (s zapasom)
+         // channel 25ps
+         qt= 50.*al*gain[i]/ph2Mip;  // 50mv/Mip amp in mV 
+         //  fill TDC
+         trCFD = Int_t (timeGaus[i] + 1000.*timeDelayCFD[i])/channelWidth; 
+         trLED= Int_t (timeGaus[i] + 1000.*timeDelayLED[i]); 
+         sl = ((TGraph*)slewingLED.At(i))->Eval(qt);
+         trLED = Int_t(( trLED + 1000*sl )/channelWidth);
+         qtCh=Int_t (1000.*TMath::Log(qt)) / channelWidth;
+         fADC0->AddAt(0,i);
+         fADC->AddAt(qtCh,i);
+         ftimeCFD->AddAt(Int_t (trCFD),i);
+         ftimeLED->AddAt(trLED,i); 
+         //      sumMult += Int_t ((al*gain[i]/ph2Mip)*50) ;
+         sumMult += Int_t (qt/sumMultCoeff)  ;
+        
+       AliDebug(10,Form("  pmt %i : time in ns %f time in channels %i   ",
+                       i, timeGaus[i],trCFD ));
+       AliDebug(10,Form(" qt in mV %f qt in ns %f qt in channels %i   ",qt, 
+                       TMath::Log(qt), qtCh));
        }
       } //pmt loop
-    for (Int_t im=0; im<3; im++)
-      { 
-       if (sumMult[im]>threshold){
-         Int_t sum=Int_t (TMath::Log(sumMult[im]*ph2mV) * mV2channel);
-         Int_t sumAmp=Int_t (TMath::Log(10*sumMult[im]*ph2mV) * mV2channel);
-         fSumMult->AddAt(sum,im);
-         fSumMult->AddAt(sumAmp,im+3);
-       }
-      }        
 
-    fSTART->AddDigit(bestRightTDC,bestLeftTDC,meanTime,timeDiff,fSumMult,
-                    ftimeTDC,fADC,ftimeTDCAmp,fADCAmp);
+    if (sumMult > threshold){
+      fSumMult =  Int_t (1000.* TMath::Log(Double_t(sumMult) / Double_t(sumMultCoeff))
+                        /channelWidth);
+      AliDebug(10,Form("summult mv %i   mult  in chammens %i in ps %i ", 
+                     sumMult, fSumMult, fSumMult*channelWidth));
+    }
+    //     if (  besttimeright<99999 || besttimeleft < 99999) {
+
+      fSTART->AddDigit(bestRightTDC,bestLeftTDC,meanTime,timeDiff,fSumMult,
+                      ftimeCFD,fADC,ftimeLED,fADC0);
+      //     } 
+     
+      AliDebug(10,Form(" Digits wrote bestRightTDC %i bestLeftTDC %i  meanTime %i  timeDiff %i fSumMult %i ", bestRightTDC,bestLeftTDC,meanTime,timeDiff,fSumMult ));
     pOutStartLoader->UnloadHits();
   } //input streams loop
   
@@ -281,33 +325,35 @@ void AliSTARTDigitizer::Exec(Option_t* /*option*/)
       pOutStartLoader->MakeTree("D");
       treeD = pOutStartLoader->TreeD();
     }
+    treeD->Reset();
+    fSTART  = (AliSTART*)outRL ->GetAliRun()->GetDetector("START");
+    // Make a branch in the tree 
     fSTART->MakeBranch("D");
-     treeD->Reset();
      treeD->Fill();
   
-  pOutStartLoader->WriteDigits("OVERWRITE");
-  
-  fSTART->ResetDigits();
-  pOutStartLoader->UnloadDigits();
-    
+     pOutStartLoader->WriteDigits("OVERWRITE");
+     
+     fSTART->ResetDigits();
+     pOutStartLoader->UnloadDigits();
+     
 }
 
 
 //------------------------------------------------------------------------
-Bool_t AliSTARTDigitizer::RegisterPhotoE(Double_t energy)
+Bool_t AliSTARTDigitizer::RegisterPhotoE(Int_t ipmt,Double_t energy)
 {
 
   
   //  Float_t hc=197.326960*1.e6; //mev*nm
   Double_t hc=1.973*1.e-6; //gev*nm
-  //  cout<<"AliSTARTDigitizer::RegisterPhotoE >> energy "<<energy<<endl;
   Float_t lambda=hc/energy;
-  Int_t bin=  fEff->GetXaxis()->FindBin(lambda);
-  Float_t eff=fEff->GetBinContent(bin);
+  Float_t eff = ((TGraph*) fEffPMT.At(ipmt))->Eval(lambda);
   Double_t  p = gRandom->Rndm();
+
   if (p > eff)
     return kFALSE;
   
   return kTRUE;
 }
+
 //----------------------------------------------------------------------------