#include <TFile.h>
//AliRoot includes
+#include "AliLog.h"
#include "AliRawReader.h"
#include "AliRawReaderRoot.h"
#include "AliRawReaderDate.h"
#include "AliRawEventHeaderBase.h"
#include "AliTPCRawStream.h"
+#include "AliTPCRawStreamFast.h"
#include "AliTPCcalibDB.h"
#include "AliTPCCalROC.h"
#include "AliTPCCalPad.h"
fNbinsRMS(100),
fXminRMS(0.1),
fXmaxRMS(5.1),
+ fPeakMinus(2),
+ fPeakPlus(3),
+ fNoiseThresholdMax(5.),
+ fNoiseThresholdSum(8.),
+ fIsZeroSuppressed(kFALSE),
fLastSector(-1),
- fOldRCUformat(kTRUE),
+ fSecRejectRatio(.4),
fROC(AliTPCROC::Instance()),
+ fMapping(NULL),
fParam(new AliTPCParam),
fPedestalTPC(0x0),
fPadNoiseTPC(0x0),
fPedestalROC(0x0),
fPadNoiseROC(0x0),
fCalRocArrayT0(72),
+ fCalRocArrayT0Err(72),
fCalRocArrayQ(72),
fCalRocArrayRMS(72),
fCalRocArrayOutliers(72),
fHistoQArray(72),
fHistoT0Array(72),
fHistoRMSArray(72),
+ fMeanT0rms(0),
+ fMeanQrms(0),
+ fMeanRMSrms(0),
fHistoTmean(72),
fParamArrayEventPol1(72),
fParamArrayEventPol2(72),
// AliTPCSignal default constructor
//
// fHTime0 = new TH1F("hTime0Event","hTime0Event",(fLastTimeBin-fFirstTimeBin)*10,fFirstTimeBin,fLastTimeBin);
+ fParam->Update();
}
//_____________________________________________________________________
AliTPCCalibCE::AliTPCCalibCE(const AliTPCCalibCE &sig) :
fNbinsRMS(sig.fNbinsRMS),
fXminRMS(sig.fXminRMS),
fXmaxRMS(sig.fXmaxRMS),
+ fPeakMinus(sig.fPeakMinus),
+ fPeakPlus(sig.fPeakPlus),
+ fNoiseThresholdMax(sig.fNoiseThresholdMax),
+ fNoiseThresholdSum(sig.fNoiseThresholdSum),
+ fIsZeroSuppressed(sig.fIsZeroSuppressed),
fLastSector(-1),
- fOldRCUformat(kTRUE),
+ fSecRejectRatio(.4),
fROC(AliTPCROC::Instance()),
+ fMapping(NULL),
fParam(new AliTPCParam),
fPedestalTPC(0x0),
fPadNoiseTPC(0x0),
fPedestalROC(0x0),
fPadNoiseROC(0x0),
fCalRocArrayT0(72),
+ fCalRocArrayT0Err(72),
fCalRocArrayQ(72),
fCalRocArrayRMS(72),
fCalRocArrayOutliers(72),
fHistoQArray(72),
fHistoT0Array(72),
fHistoRMSArray(72),
+ fMeanT0rms(sig.fMeanT0rms),
+ fMeanQrms(sig.fMeanQrms),
+ fMeanRMSrms(sig.fMeanRMSrms),
fHistoTmean(72),
fParamArrayEventPol1(72),
fParamArrayEventPol2(72),
fVEventTime.SetElements(sig.fVEventTime.GetMatrixArray());
fVEventNumber.SetElements(sig.fVEventNumber.GetMatrixArray());
+ fParam->Update();
}
//_____________________________________________________________________
AliTPCCalibCE& AliTPCCalibCE::operator = (const AliTPCCalibCE &source)
//
fCalRocArrayT0.Delete();
+ fCalRocArrayT0Err.Delete();
fCalRocArrayQ.Delete();
fCalRocArrayRMS.Delete();
fCalRocArrayOutliers.Delete();
if ( fDebugStreamer) delete fDebugStreamer;
// if ( fHTime0 ) delete fHTime0;
- delete fROC;
+// delete fROC;
delete fParam;
}
//_____________________________________________________________________
// no extra analysis necessary. Assumes knowledge of the signal shape!
// assumes that it is looped over consecutive time bins of one pad
//
+
+ //temp
+// if (icsector<36) return 0;
+// if (icsector%36>17) return 0;
+
+
+ if (icRow<0) return 0;
+ if (icPad<0) return 0;
+ if (icTimeBin<0) return 0;
if ( (icTimeBin>fLastTimeBin) || (icTimeBin<fFirstTimeBin) ) return 0;
Int_t iChannel = fROC->GetRowIndexes(icsector)[icRow]+icPad; // global pad position in sector
}
if ( fPedestalROC&&fPadNoiseROC ){
- fPadPedestal = fPedestalROC->GetValue(fCurrentChannel);
+ fPadPedestal = fPedestalROC->GetValue(fCurrentChannel)*(Float_t)(!fIsZeroSuppressed);
fPadNoise = fPadNoiseROC->GetValue(fCurrentChannel);
noPedestal = kFALSE;
}
//if we are not running with pedestal database, or for the current sector there is no information
//available, calculate the pedestal and noise on the fly
if ( noPedestal ) {
+ fPadPedestal = 0;
+ fPadNoise = 0;
+ if ( fIsZeroSuppressed ) return;
const Int_t kPedMax = 100; //maximum pedestal value
Float_t max = 0;
Float_t maxPos = 0;
rms +=histo[median+idelta]*(median+idelta)*(median+idelta);
}
}
- fPadPedestal = 0;
- fPadNoise = 0;
if ( count > 0 ) {
mean/=count;
rms = TMath::Sqrt(TMath::Abs(rms/count-mean*mean));
Float_t ceQmax =0, ceQsum=0, ceTime=0, ceRMS=0;
Int_t cemaxpos = 0;
- Float_t ceSumThreshold = 8.*fPadNoise; // threshold for the signal sum
+ Float_t ceSumThreshold = fNoiseThresholdSum*fPadNoise; // threshold for the signal sum
const Int_t kCemin = 4; // range for the analysis of the ce signal +- channels from the peak
const Int_t kCemax = 7;
//
// Find local maxima on the pad signal and Histogram them
//
- Float_t ceThreshold = 5.*fPadNoise; // threshold for the signal
+ Float_t ceThreshold = fNoiseThresholdMax*TMath::Max(fPadNoise,Float_t(1.)); // threshold for the signal
Int_t count = 0;
- Int_t tminus = 2;
- Int_t tplus = 3;
- for (Int_t i=fLastTimeBin-tplus-1; i>=fFirstTimeBin+tminus; --i){
- if ( (fPadSignal[i]-fPadPedestal)>ceThreshold && IsPeak(i,tminus,tplus) ){
+// Int_t tminus = 2;
+// Int_t tplus = 3;
+ for (Int_t i=fLastTimeBin-fPeakPlus-1; i>=fFirstTimeBin+fPeakMinus; --i){
+ if ( (fPadSignal[i]-fPadPedestal)>ceThreshold && IsPeak(i,fPeakMinus,fPeakPlus) ){
if (count<maxima.GetNrows()){
maxima.GetMatrixArray()[count++]=i;
GetHistoTmean(fCurrentSector,kTRUE)->Fill(i);
FindLocalMaxima(maxima);
if ( (fNevents == 0) || (fOldRunNumber!=fRunNumber) ) return; // return because we don't have Time0 info for the CE yet
-
+ if ( !GetTMeanEvents(fCurrentSector) ) return; //return if we don't have time 0 info, eg if only one side has laser
TVectorD param(3);
Float_t qSum;
//check if last pad has allready been processed, if not do so
if ( fMaxTimeBin>-1 ) ProcessPad();
+ AliDebug(3, Form("EndEvent() - Start; Event: %05d", fNevents));
+
TVectorD param(3);
TMatrixD dummy(3,3);
// TVectorF vMeanTime(72);
AliTPCCalROC *calOroc=new AliTPCCalROC(36);
//find mean time0 offset for side A and C
- Double_t time0Side[2]; //time0 for side A:0 and C:0
- Double_t time0SideCount[2]; //time0 counter for side A:0 and C:0
+ Double_t time0Side[2]; //time0 for side A:0 and C:1
+ Double_t time0SideCount[2]; //time0 counter for side A:0 and C:1
time0Side[0]=0;time0Side[1]=0;time0SideCount[0]=0;time0SideCount[1]=0;
for ( Int_t iSec = 0; iSec<72; ++iSec ){
time0Side[(iSec/18)%2] += fVTime0Offset.GetMatrixArray()[iSec];
time0Side[1]/=time0SideCount[1];
// end find time0 offset
+ Int_t nSecMeanT=0;
//loop over all ROCs, fill CE Time histogram corrected for the mean Time0 of each ROC
for ( Int_t iSec = 0; iSec<72; ++iSec ){
+ AliDebug(4,Form("Processing sector '%02d'\n",iSec));
//find median and then calculate the mean around it
TH1S *hMeanT = GetHistoTmean(iSec); //histogram with local maxima position information
if ( !hMeanT ) continue;
+ //continue if not enough data is filled in the meanT histogram. This is the case if we do not have a laser event.
+ if ( hMeanT->GetEntries() < fROC->GetNChannels(iSec)*fSecRejectRatio ){
+ hMeanT->Reset();
+ AliDebug(3,Form("Skipping sec. '%02d': Not enough statistics\n",iSec));
+ continue;
+ }
Double_t entries = hMeanT->GetEntries();
Double_t sum = 0;
Int_t ibin=0;
for ( ibin=0; ibin<hMeanT->GetNbinsX(); ++ibin){
sum+=arr[ibin];
- if ( sum>=(entries/2) ) break;
+ if ( sum>=(entries/2.) ) break;
}
Int_t delta = 4;
Int_t firstBin = fFirstTimeBin+ibin-delta;
vMeanTime->ResizeTo(vSize+100);
vMeanTime->GetMatrixArray()[fNevents]=median;
+ nSecMeanT++;
// end find median
-
+
TVectorF *vTimes = GetPadTimesEvent(iSec);
if ( !vTimes ) continue; //continue if no time information for this sector is available
}
//----------------------------- Debug end ------------------------------
}// end channel loop
- hMeanT->Reset();
TVectorD paramPol1(3);
TVectorD paramPol2(6);
GetParamArrayPol1(iSec,kTRUE)->AddAtAndExpand(new TVectorD(paramPol1), fNevents);
GetParamArrayPol2(iSec,kTRUE)->AddAtAndExpand(new TVectorD(paramPol2), fNevents);
}
-// printf("events: %d -- size: %d\n",fNevents,GetParamArrayPol1(iSec)->GetSize());
//------------------------------- Debug start ------------------------------
if ( fDebugLevel>0 ){
"\n";
}
//------------------------------- Debug end ------------------------------
+ hMeanT->Reset();
}// end sector loop
-
+ //return if no sector has a valid mean time
+ if ( nSecMeanT == 0 ) return;
+
+
// fTMeanArrayEvent.AddAtAndExpand(new TVectorF(vMeanTime),fNevents);
// fQMeanArrayEvent.AddAtAndExpand(new TVectorF(vMeanQ),fNevents);
if ( fVEventTime.GetNrows() < fNevents+1 ) {
delete calIroc;
delete calOroc;
+ AliDebug(3, Form("EndEvent() - End; Event: %05d", fNevents));
+}
+//_____________________________________________________________________
+Bool_t AliTPCCalibCE::ProcessEventFast(AliTPCRawStreamFast *rawStreamFast)
+{
+ //
+ // Event Processing loop - AliTPCRawStreamFast
+ //
+ ResetEvent();
+ Bool_t withInput = kFALSE;
+ while ( rawStreamFast->NextDDL() ){
+ while ( rawStreamFast->NextChannel() ){
+ Int_t isector = rawStreamFast->GetSector(); // current sector
+ Int_t iRow = rawStreamFast->GetRow(); // current row
+ Int_t iPad = rawStreamFast->GetPad(); // current pad
+
+ while ( rawStreamFast->NextBunch() ){
+ Int_t startTbin = (Int_t)rawStreamFast->GetStartTimeBin();
+ Int_t endTbin = (Int_t)rawStreamFast->GetEndTimeBin();
+ for (Int_t iTimeBin = startTbin; iTimeBin < endTbin; iTimeBin++){
+ Float_t signal=(Float_t)rawStreamFast->GetSignals()[iTimeBin-startTbin];
+ Update(isector,iRow,iPad,iTimeBin+1,signal);
+ withInput = kTRUE;
+ }
+ }
+ }
+ }
+ if (withInput){
+ EndEvent();
+ }
+ return withInput;
+}
+//_____________________________________________________________________
+Bool_t AliTPCCalibCE::ProcessEventFast(AliRawReader *rawReader)
+{
+ //
+ // Event processing loop using the fast raw stream algorithm- AliRawReader
+ //
+
+ //printf("ProcessEventFast - raw reader\n");
+
+ AliRawEventHeaderBase* eventHeader = (AliRawEventHeaderBase*)rawReader->GetEventHeader();
+ if (eventHeader){
+ fTimeStamp = eventHeader->Get("Timestamp");
+ fRunNumber = eventHeader->Get("RunNb");
+ }
+ fEventId = *rawReader->GetEventId();
+
+ AliTPCRawStreamFast *rawStreamFast = new AliTPCRawStreamFast(rawReader, (AliAltroMapping**)fMapping);
+ Bool_t res=ProcessEventFast(rawStreamFast);
+ delete rawStreamFast;
+ return res;
+
}
//_____________________________________________________________________
Bool_t AliTPCCalibCE::ProcessEvent(AliTPCRawStream *rawStream)
// The Function 'SetTimeStamp' should be called for each event to set the event time stamp!!!
//
- rawStream->SetOldRCUFormat(fOldRCUformat);
-
ResetEvent();
Bool_t withInput = kFALSE;
//
- AliTPCRawStream rawStream(rawReader);
+ AliTPCRawStream rawStream(rawReader,(AliAltroMapping**)fMapping);
AliRawEventHeaderBase* eventHeader = (AliRawEventHeaderBase*)rawReader->GetEventHeader();
if (eventHeader){
fTimeStamp = eventHeader->Get("Timestamp");
sprintf(name,"hCalib%s%.2d",type,sector);
sprintf(title,"%s calibration histogram sector %.2d",type,sector);
- // new histogram with Q calib information. One value for each pad!
+ // new histogram with calib information. One value for each pad!
TH1S* hist = new TH1S(name,title,
fLastTimeBin-fFirstTimeBin,fFirstTimeBin,fLastTimeBin);
hist->SetDirectory(0);
AliTPCCalROC* AliTPCCalibCE::GetCalRocT0(Int_t sector, Bool_t force)
{
//
- // return pointer to Carge ROC Calibration
+ // return pointer to Time 0 ROC Calibration
// if force is true create a new histogram if it doesn't exist allready
//
TObjArray *arr = &fCalRocArrayT0;
return GetCalRoc(sector, arr, force);
}
//_____________________________________________________________________
+AliTPCCalROC* AliTPCCalibCE::GetCalRocT0Err(Int_t sector, Bool_t force)
+{
+ //
+ // return pointer to the error of Time 0 ROC Calibration
+ // if force is true create a new histogram if it doesn't exist allready
+ //
+ TObjArray *arr = &fCalRocArrayT0Err;
+ return GetCalRoc(sector, arr, force);
+}
+//_____________________________________________________________________
AliTPCCalROC* AliTPCCalibCE::GetCalRocQ(Int_t sector, Bool_t force)
{
//
arrPol2->Expand(fNevents+nCEevents);
}
if ( vMeanTimeCE && vMeanQCE ){
- vMeanTime = GetTMeanEvents(iSec);
- vMeanQCE = GetQMeanEvents(iSec);
+ vMeanTime = GetTMeanEvents(iSec,kTRUE);
+ vMeanQ = GetQMeanEvents(iSec,kTRUE);
vMeanTime->ResizeTo(fNevents+nCEevents);
vMeanQ->ResizeTo(fNevents+nCEevents);
}
TVectorD paramRMS(3);
TMatrixD dummy(3,3);
+ Float_t channelCounter=0;
+ fMeanT0rms=0;
+ fMeanQrms=0;
+ fMeanRMSrms=0;
+
for (Int_t iSec=0; iSec<72; ++iSec){
TH2S *hT0 = GetHistoT0(iSec);
if (!hT0 ) continue;
- AliTPCCalROC *rocQ = GetCalRocQ (iSec,kTRUE);
- AliTPCCalROC *rocT0 = GetCalRocT0 (iSec,kTRUE);
- AliTPCCalROC *rocRMS = GetCalRocRMS(iSec,kTRUE);
- AliTPCCalROC *rocOut = GetCalRocOutliers(iSec,kTRUE);
+ AliTPCCalROC *rocQ = GetCalRocQ (iSec,kTRUE);
+ AliTPCCalROC *rocT0 = GetCalRocT0 (iSec,kTRUE);
+ AliTPCCalROC *rocT0Err = GetCalRocT0Err (iSec,kTRUE);
+ AliTPCCalROC *rocRMS = GetCalRocRMS(iSec,kTRUE);
+ AliTPCCalROC *rocOut = GetCalRocOutliers(iSec,kTRUE);
TH2S *hQ = GetHistoQ(iSec);
TH2S *hRMS = GetHistoRMS(iSec);
Float_t cogTime0 = -1000;
Float_t cogQ = -1000;
Float_t cogRMS = -1000;
- Float_t cogOut = 0;
+ Float_t cogOut = 0;
+ Float_t rms = 0;
+ Float_t rmsT0 = 0;
Int_t offsetQ = (fNbinsQ+2)*(iChannel+1)+1;
Int_t offsetT0 = (fNbinsT0+2)*(iChannel+1)+1;
Int_t offsetRMS = (fNbinsRMS+2)*(iChannel+1)+1;
- cogQ = AliMathBase::GetCOG(arrayhQ+offsetQ,fNbinsQ,fXminQ,fXmaxQ);
- cogTime0 = AliMathBase::GetCOG(arrayhT0+offsetT0,fNbinsT0,fXminT0,fXmaxT0);
- cogRMS = AliMathBase::GetCOG(arrayhRMS+offsetRMS,fNbinsRMS,fXminRMS,fXmaxRMS);
-
-
+ cogQ = AliMathBase::GetCOG(arrayhQ+offsetQ,fNbinsQ,fXminQ,fXmaxQ,&rms);
+ fMeanQrms+=rms;
+ cogTime0 = AliMathBase::GetCOG(arrayhT0+offsetT0,fNbinsT0,fXminT0,fXmaxT0,&rmsT0);
+ fMeanT0rms+=rmsT0;
+ cogRMS = AliMathBase::GetCOG(arrayhRMS+offsetRMS,fNbinsRMS,fXminRMS,fXmaxRMS,&rms);
+ fMeanRMSrms+=rms;
+ channelCounter++;
/*
//outlier specifications
*/
rocQ->SetValue(iChannel, cogQ*cogQ);
rocT0->SetValue(iChannel, cogTime0);
+ rocT0Err->SetValue(iChannel, rmsT0);
rocRMS->SetValue(iChannel, cogRMS);
rocOut->SetValue(iChannel, cogOut);
}
}
+ if ( channelCounter>0 ){
+ fMeanT0rms/=channelCounter;
+ fMeanQrms/=channelCounter;
+ fMeanRMSrms/=channelCounter;
+ }
if ( fDebugStreamer ) fDebugStreamer->GetFile()->Write();
// delete fDebugStreamer;
// fDebugStreamer = 0x0;
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff