1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
18 ////////////////////////////////////////////////////////////////////////////////////////
20 // Implementation of the TPC Raw drift velocity and Altro L1 Phase calibration //
22 // Origin: Jens Wiechula, J.Wiechula@gsi.de //
24 ////////////////////////////////////////////////////////////////////////////////////////
27 // *************************************************************************************
28 // * Class Description *
29 // *************************************************************************************
33 TFile f("CalibAltro.root");
34 AliTPCCalibRaw *al=(AliTPCCalibRaw*)f.Get(f.GetListOfKeys()->At(0)->GetName())
36 TCanvas *c1=(TCanvas*)gROOT->FindObject("c1");
37 if (!c1) c1=new TCanvas("c1","c1");
40 TH2F h2f("h2","h2;RCU;fail",216,0,216,al->GetNevents(),0,al->GetNevents());
44 TVectorF *v=al->GetALTROL1PhaseFailEventsRCU(i);
46 for (iev=0;iev<al->GetNevents();++iev) {
47 h2f->SetBinContent(i+1,iev+1,(*v)(iev));
50 // h.SetLineColor(i/216.*50+50);
51 // ((TH1F*)h.Clone(Form("h%d",i)))->Draw(first?"":"same");
69 #include <TObjString.h>
70 #include <TTimeStamp.h>
73 #include "AliTPCCalROC.h"
74 #include "AliAltroRawStream.h"
77 #include "AliTPCCalibRaw.h"
79 ClassImp(AliTPCCalibRaw)
81 AliTPCCalibRaw::AliTPCCalibRaw() :
86 fNFailL1PhaseEvent(0),
88 fNSecTime(600), //default 10 minutes
89 fNBinsTime(60), //default 60*10 minutes = 10 hours
90 fPadProcessed(kFALSE),
102 fArrCurrentPhaseDist(4),
103 fArrCurrentPhase(kNRCU),
104 fArrFailEventNumber(100),
105 fArrALTROL1Phase(100000),
106 fArrALTROL1PhaseEvent(kNRCU),
107 fArrALTROL1PhaseFailEvent(kNRCU),
109 fVOccupancyEvent(100000),
110 fVSignalSumEvent(100000),
111 fVOccupancySenEvent(100000),
112 fVSignalSumSenEvent(100000),
113 fVNfiredPadsSenEvent(100000),
114 fVTimeStampEvent(100000)
119 SetNameTitle("AliTPCCalibRaw","AliTPCCalibRaw");
121 for (Int_t ircu=0;ircu<kNRCU;++ircu) fArrCurrentPhase.GetMatrixArray()[ircu]=-1;
125 //_____________________________________________________________________
126 AliTPCCalibRaw::AliTPCCalibRaw(const TMap *config) :
127 AliTPCCalibRawBase(),
131 fNFailL1PhaseEvent(0),
133 fNSecTime(600), //default 10 minutes
134 fNBinsTime(60), //default 60*10 minutes = 10 hours
135 fPadProcessed(kFALSE),
147 fArrCurrentPhaseDist(4),
148 fArrCurrentPhase(kNRCU),
149 fArrFailEventNumber(100),
150 fArrALTROL1Phase(100000),
151 fArrALTROL1PhaseEvent(kNRCU),
152 fArrALTROL1PhaseFailEvent(kNRCU),
154 fVOccupancyEvent(100000),
155 fVSignalSumEvent(100000),
156 fVOccupancySenEvent(100000),
157 fVSignalSumSenEvent(100000),
158 fVNfiredPadsSenEvent(100000),
159 fVTimeStampEvent(100000)
164 SetNameTitle("AliTPCCalibRaw","AliTPCCalibRaw");
166 for (Int_t ircu=0;ircu<kNRCU;++ircu) fArrCurrentPhase.GetMatrixArray()[ircu]=-1;
169 if (config->GetValue("FirstTimeBin")) fFirstTimeBin = ((TObjString*)config->GetValue("FirstTimeBin"))->GetString().Atoi();
170 if (config->GetValue("LastTimeBin")) fLastTimeBin = ((TObjString*)config->GetValue("LastTimeBin"))->GetString().Atoi();
171 if (config->GetValue("DebugLevel")) fDebugLevel = ((TObjString*)config->GetValue("DebugLevel"))->GetString().Atoi();
174 //_____________________________________________________________________
175 AliTPCCalibRaw::~AliTPCCalibRaw()
182 //_____________________________________________________________________
183 // AliTPCCalibRaw& AliTPCCalibRaw::operator = (const AliTPCCalibRaw &source)
186 // // assignment operator
188 // if (&source == this) return *this;
189 // new (this) AliTPCCalibRaw(source);
194 //_____________________________________________________________________
195 Int_t AliTPCCalibRaw::Update(const Int_t isector, const Int_t iRow, const Int_t iPad,
196 const Int_t iTimeBin, const Float_t signal)
199 // Data filling method
201 if (iRow<0) return 0;
202 if (iPad<0) return 0;
203 if (iTimeBin<0) return 0;
204 if (!fFirstTimeStamp) fFirstTimeStamp=GetTimeStamp();
206 Int_t iChannel = fROC->GetRowIndexes(isector)[iRow]+iPad; // global pad position in sector
208 fVOccupancyEvent.GetMatrixArray()[GetNevents()]++;
209 fVSignalSumEvent.GetMatrixArray()[GetNevents()]+=signal;
210 //occupancy in sensitive regions
211 Int_t npads=(Int_t)fROC->GetNPads(isector,iRow);
212 Int_t cpad=iPad-npads/2;
213 if (isector<(Int_t)fROC->GetNInnerSector()){
215 if ( iRow>19 && iRow<46 ){
216 if ( TMath::Abs(cpad)<7 ){
217 fVOccupancySenEvent.GetMatrixArray()[GetNevents()]++;
218 fVSignalSumSenEvent.GetMatrixArray()[GetNevents()]+=signal;
219 if (iChannel!=fCurrentChannel) fVNfiredPadsSenEvent.GetMatrixArray()[GetNevents()]++;
222 } else if ( iRow>75 ){
223 //OROC case (outer corners and last three rows are sensitive)
224 Int_t padEdge=(Int_t)TMath::Min(iPad,npads-iPad);
225 Int_t nrows=(Int_t)fROC->GetNRows(isector);
226 if ((nrows-iRow-1)<3 || padEdge<((((Int_t)iRow-76)/4+1))*2){
227 fVOccupancySenEvent.GetMatrixArray()[GetNevents()]++;
228 fVSignalSumSenEvent.GetMatrixArray()[GetNevents()]+=signal;
229 if (iChannel!=fCurrentChannel) fVNfiredPadsSenEvent.GetMatrixArray()[GetNevents()]++;
233 if ( (iTimeBin>fLastTimeBin) || (iTimeBin<fFirstTimeBin) ) return 0;
234 //don't process edge pads
235 if (IsEdgePad(isector,iRow,iPad)) return 0;
236 // Double_t x[kHnBinsDV]={1,isector,0};
237 // fHnDrift->Fill(x);
238 if (fCurrentChannel==iChannel){
239 if (fPadProcessed) return 0;
241 fPadProcessed=kFALSE;
247 // Double_t x2[kHnBinsDV]={2,isector,0};
248 // fHnDrift->Fill(x2);
251 if (signal>fLastSignal) ++fNOkPlus;
252 else if(signal<fLastSignal && fNOkPlus>=fPeakDetPlus){
254 if (!fPeakTimeBin) fPeakTimeBin=fLastTimeBinProc;
255 if ( fNOkMinus>=fPeakDetMinus ) {
256 Double_t x[kHnBinsDV]={fPeakTimeBin,isector,(fTimeStamp-fFirstTimeStamp)/fNSecTime};
266 fLastTimeBinProc=iTimeBin;
267 fLastSignal=TMath::Nint(signal);
268 fCurrentChannel = iChannel;
271 //_____________________________________________________________________
272 void AliTPCCalibRaw::UpdateDDL(){
274 // fill ALTRO L1 information
280 fNanoSec=s.GetNanoSec();
283 Int_t phase=(Int_t)(GetL1PhaseTB()*4.);
284 //Fill pahse information of current rcu and event
285 fArrCurrentPhase.GetMatrixArray()[fCurrDDLNum]=phase;
286 //increase phase counter
287 ++((fArrCurrentPhaseDist.GetMatrixArray())[phase]);
290 //_____________________________________________________________________
291 void AliTPCCalibRaw::ResetEvent()
294 // Reset event counters
300 fArrCurrentPhaseDist.Zero();
302 //_____________________________________________________________________
303 void AliTPCCalibRaw::EndEvent()
306 // End event analysis
310 //find phase of the current event
311 Int_t phaseMaxEntries=-1;
313 for (Int_t i=0;i<fArrCurrentPhaseDist.GetNrows();++i){
314 Int_t entries=(Int_t)fArrCurrentPhaseDist[i];
315 if (maxEntries<entries) {
320 // store phase of current event
321 if (fArrALTROL1Phase.GetNrows()-1<=GetNevents())
322 fArrALTROL1Phase.ResizeTo(GetNevents()+10000);
323 (fArrALTROL1Phase.GetMatrixArray())[GetNevents()]=phaseMaxEntries;
325 //loop over RCUs and test failures
327 for (Int_t ircu=0;ircu<kNRCU;++ircu){
328 Int_t phase=(Int_t)fArrCurrentPhase[ircu];
329 if (phase<0) continue;
330 if (phase!=phaseMaxEntries){
331 TVectorF *arr=MakeArrL1PhaseRCU(fCurrDDLNum,kTRUE);
332 if (arr->GetNrows()-1<=(Int_t)fNFailL1PhaseEvent) arr->ResizeTo(arr->GetNrows()+100);
333 (arr->GetMatrixArray())[fNFailL1PhaseEvent]=phase;
337 //reset current phase information
338 fArrCurrentPhase[ircu]=-1;
341 if (fArrFailEventNumber.GetNrows()-1<=(Int_t)fNFailL1PhaseEvent) fArrFailEventNumber.ResizeTo(fArrFailEventNumber.GetNrows()+100);
342 fArrFailEventNumber.GetMatrixArray()[fNFailL1PhaseEvent]=GetNevents();
344 fNFailL1PhaseEvent+=fail;
346 fVTimeStampEvent.GetMatrixArray()[GetNevents()]=GetTimeStamp()-fFirstTimeStamp+fNanoSec*1e-9;
349 if (fVOccupancyEvent.GetNrows()-1<=GetNevents()){
350 fVOccupancyEvent.ResizeTo(GetNevents()+10000);
351 fVSignalSumEvent.ResizeTo(GetNevents()+10000);
352 fVOccupancySenEvent.ResizeTo(GetNevents()+10000);
353 fVSignalSumSenEvent.ResizeTo(GetNevents()+10000);
354 fVTimeStampEvent.ResizeTo(GetNevents()+10000);
355 fVNfiredPadsSenEvent.ResizeTo(GetNevents()+10000);
359 //_____________________________________________________________________
360 TH2C *AliTPCCalibRaw::MakeHistL1RCUEvents(Int_t type)
362 // Create a 2D histo RCU:Events indicating the there was a deviation
363 // from the mean L1 phase of the event
365 //type: 0=Failures, 1=Phases
367 //number of relavant events, depending on version
368 Int_t nevents=GetNevents();
370 Bool_t newVersion=kFALSE;
371 for (Int_t ircu=0; ircu<kNRCU; ++ircu){
372 const TVectorF *v=GetALTROL1PhaseEventsRCU(ircu);
374 if ((UInt_t)(v->GetNrows())==fNFailL1PhaseEvent){
376 nevents=fNFailL1PhaseEvent;
380 TH2C *h2 = new TH2C("hL1FailRCUEvents","L1 Failures;RCU;Event",kNRCU,0,kNRCU,nevents,0,nevents);
382 for (Int_t ircu=0;ircu<kNRCU;++ircu) {
383 const TVectorF *v=GetALTROL1PhaseEventsRCU(ircu);
388 } else if (type==1) {
394 for (Int_t iev=0;iev<nevents;++iev) {
395 Float_t val=(*v)(iev);
396 Float_t phase=fArrALTROL1Phase.GetMatrixArray()[iev];
398 Int_t event=(Int_t)fArrFailEventNumber.GetMatrixArray()[iev];
399 phase=fArrALTROL1Phase.GetMatrixArray()[event];
401 if (type==0) val=(val!=phase);
402 h2->SetBinContent(ircu+1,iev+1,val+add);
407 //_____________________________________________________________________
408 TH1F *AliTPCCalibRaw::MakeHistL1PhaseDist()
411 // L1 phase distribution. Should be flat in ideal case
413 TH1F *h=new TH1F("L1phaseDist","Normalized L1 phase distribution;phase;fraction of events",4,0,4);
415 for (Int_t iev=0;iev<GetNevents();++iev) h->Fill(fArrALTROL1Phase.GetMatrixArray()[iev]);
416 if (GetNevents()>0) h->Scale(1./GetNevents());
421 //_____________________________________________________________________
422 TVectorF *AliTPCCalibRaw::MakeVectL1PhaseDist()
425 // L1 phase distribution. Should be flat in ideal case
427 TVectorF *v=new TVectorF(4);
428 for (Int_t iev=0;iev<GetNevents();++iev) {
429 Int_t phase=(Int_t)fArrALTROL1Phase.GetMatrixArray()[iev];
430 ((v->GetMatrixArray())[phase])+=1./GetNevents();
434 //_____________________________________________________________________
435 TH2C *AliTPCCalibRaw::MakeHistL1RCUEventsIROC(Int_t type)
438 // Create a 2D histo RCU:Events indicating the there was a deviation
439 // from the mean L1 phase of the event
441 TH2C *h2 = new TH2C("hL1FailRCUEventsIROC","L1 Failures IROCs;RCU;Event",72,0,36,GetNevents(),0,GetNevents());
442 for (Int_t ircu=0;ircu<72;++ircu) {
444 if (type==0) v=GetALTROL1PhaseFailEventsRCU(ircu);
445 else if (type==1) v=GetALTROL1PhaseEventsRCU(ircu);
447 for (Int_t iev=0;iev<GetNevents();++iev) {
448 h2->SetBinContent(ircu+1,iev+1,(*v)(iev));
453 //_____________________________________________________________________
454 TH2C *AliTPCCalibRaw::MakeHistL1RCUEventsOROC(Int_t type)
457 // Create a 2D histo RCU:Events indicating the there was a deviation
458 // from the mean L1 phase of the event
460 TH2C *h2 = new TH2C("hL1FailRCUEventsOROC","L1 Failures OROCs;RCU;Event",144,0,36,GetNevents(),0,GetNevents());
461 for (Int_t ircu=72;ircu<kNRCU;++ircu) {
463 if (type==0) v=GetALTROL1PhaseFailEventsRCU(ircu);
464 else if (type==1) v=GetALTROL1PhaseEventsRCU(ircu);
466 for (Int_t iev=0;iev<GetNevents();++iev) {
467 h2->SetBinContent(ircu-72+1,iev+1,(*v)(iev));
472 //_____________________________________________________________________
473 void AliTPCCalibRaw::CreateDVhist()
476 // Setup the HnSparse for the drift velocity determination
478 if (fHnDrift) return;
480 //time bin, roc, time
481 Int_t bins[kHnBinsDV] = {fLastTimeBin-fFirstTimeBin, 72, fNBinsTime};
482 Double_t xmin[kHnBinsDV] = {fFirstTimeBin,0,0};
483 Double_t xmax[kHnBinsDV] = {fLastTimeBin,72,fNBinsTime};
484 fHnDrift=new THnSparseI("fHnDrift",Form("Drift velocity using last time bin;time bin[#times 100ns];ROC;Time bin [#times %us]",fNSecTime),kHnBinsDV, bins, xmin, xmax);
487 //_____________________________________________________________________
488 void AliTPCCalibRaw::Analyse()
495 fArrALTROL1Phase.ResizeTo(GetNevents());
496 for (Int_t ircu=0;ircu<kNRCU;++ircu){
497 TVectorF *arr=MakeArrL1PhaseRCU(ircu);//MakeArrL1PhaseRCU(ircu);
499 arr->ResizeTo(fNFailL1PhaseEvent);
500 fArrFailEventNumber.ResizeTo(fNFailL1PhaseEvent);
501 // TVectorF *arrF=MakeArrL1PhaseFailRCU(ircu);
502 // arrF->ResizeTo(1);
504 //resize occupancy arrays only save event occupancy in sensitive regions by default
505 //save the rest in debub mode
506 fVOccupancySenEvent.ResizeTo(GetNevents());
508 fVOccupancyEvent.ResizeTo(GetNevents());
509 fVSignalSumEvent.ResizeTo(GetNevents());
510 fVSignalSumSenEvent.ResizeTo(GetNevents());
511 fVNfiredPadsSenEvent.ResizeTo(GetNevents());
512 fVTimeStampEvent.ResizeTo(GetNevents());
514 fVOccupancyEvent.ResizeTo(0);
515 fVSignalSumEvent.ResizeTo(0);
516 fVSignalSumSenEvent.ResizeTo(0);
517 fVNfiredPadsSenEvent.ResizeTo(0);
518 fVTimeStampEvent.ResizeTo(0);
520 //Analyse drift velocity TODO
523 //_____________________________________________________________________
524 TGraph* AliTPCCalibRaw::MakeGraphOccupancy(const Int_t type, const Int_t xType)
527 // create occupancy graph (of samples abouve threshold)
528 // type=0: number of samples
529 // type=1: mean data volume (ADC counts/sample)
530 // type=2: data volume (ADC counts)
531 // type=3: samples per ADC count
532 // type=4: sample occupancy
534 // type=5: number of sample sensitive / number of samples
536 // same in sensitive regions:
537 // type=10: number of samples
538 // type=11: mean data volume (ADC counts/sample)
539 // type=12: data volume (ADC counts)
540 // type=13: samples per ADC count
541 // type=14: sample occupancy
543 // type=16: number of samples sensitive / number of pads sensitive
544 // type=17: pad occupancy in sensitive regions
545 // xType=0: vs. time stamp
546 // xType=1: vs. event counter
549 TString title("Event occupancy");
550 TString xTitle("Time");
551 TString yTitle("number of samples");
552 TGraph *gr=new TGraph(GetNevents());
553 if (fVSignalSumEvent.GetNrows()==0&&!(type==10||type==14)) return 0;
554 TVectorF *vOcc=&fVOccupancyEvent;
555 TVectorF *vSum=&fVSignalSumEvent;
556 TVectorF *vPads=&fVNfiredPadsSenEvent;
557 Double_t norm=557568.;
558 if (type!=14&&fVOccupancyEvent.GetNrows()==0){
559 AliWarning("In non debug mode only occupancy in sensitive regions vs. event awailable!!!");
563 vOcc=&fVOccupancySenEvent;
564 vSum=&fVSignalSumSenEvent;
565 vPads=&fVNfiredPadsSenEvent;
568 for (Int_t i=0;i<GetNevents(); ++i){
569 Double_t nAboveThreshold=vOcc->GetMatrixArray()[i];
572 Double_t timestamp =1;
575 if (fVOccupancyEvent.GetNrows()>0){
576 nSumADC =vSum->GetMatrixArray()[i];
577 timestamp =fVTimeStampEvent.GetMatrixArray()[i]+fFirstTimeStamp;
578 nPads =vPads->GetMatrixArray()[i];
580 Double_t x=timestamp;
585 if (type%10==0) y=nAboveThreshold;
586 if (type%10==1&&nAboveThreshold>0) y=nSumADC/nAboveThreshold;
587 if (type%10==2) y=nSumADC;
588 if (type%10==3&&nSumADC>0) y=nAboveThreshold/nSumADC;
589 if (type%10==4) y=nAboveThreshold/(norm*(fLastTimeBin-fFirstTimeBin));
590 if (type==5) y=fVOccupancySenEvent.GetMatrixArray()[i]/fVOccupancyEvent.GetMatrixArray()[i];
591 if (type==16&&nPads>0) y=nAboveThreshold/nPads;
592 if (type==17) y=nPads/norm;
596 if (xType==1) xTitle="Event";
597 if (type%10==1) yTitle="Mean ADC counts/sample";
598 else if (type%10==2) yTitle="Data volume [ADC counts]";
599 else if (type%10==3) yTitle="samples per ADC count";
600 else if (type%10==4) yTitle="sample occupancy";
601 if (type==5) yTitle="N samples (sensitive) / N samples";
602 if (type%10==6) yTitle="N samples / N pads";
603 if (type==17) yTitle="Pad Occupancy";
604 if (type>=10) yTitle+=" (sensitive)";
605 title=yTitle+":"+xTitle;
606 title+=";"+xTitle+";"+yTitle;
607 gr->SetTitle(title.Data());
608 gr->SetEditable(kFALSE);
611 //_____________________________________________________________________
612 TGraph* AliTPCCalibRaw::MakeGraphNoiseEvents()
619 //_____________________________________________________________________
620 TCanvas* AliTPCCalibRaw::MakeCanvasOccupancy(const Int_t xType, Bool_t sen)
623 // Create a canvas with occupancy information of all 'type's (see MakeGraphOccupancy)
624 // xType=0: vs. timestamp
625 // xType=1: vs. event number
627 // sen=kTRUE: for sensitive regions
630 TString name("RawOccupancy_");
631 TString title("Raw Occupancy vs. ");
635 } else if (xType==1){
641 title+=" (sensitive)";
643 TCanvas *c=(TCanvas*)gROOT->GetListOfCanvases()->FindObject(name.Data());
644 if (!c) c=new TCanvas(name.Data(),title.Data());
647 for (Int_t i=0;i<4;++i){
649 TGraph *gr=MakeGraphOccupancy(i+10*(Int_t)sen,xType);