1 /**************************************************************************
4 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
10 * Author: The ALICE Off-line Project. *
13 * Contributors are mentioned in the code where appropriate. *
19 * Permission to use, copy, modify and distribute this software and its *
22 * documentation strictly for non-commercial purposes is hereby granted *
25 * without fee, provided that the above copyright notice appears in all *
28 * copies and that both the copyright notice and this permission notice *
31 * appear in the supporting documentation. The authors make no claims *
34 * about the suitability of this software for any purpose. It is *
37 * provided "as is" without express or implied warranty. *
40 **************************************************************************/
46 /* $Id: AliTPCclusterKr.cxx,v 1.7 2008/01/22 17:24:53 matyja Exp $ */
52 //-----------------------------------------------------------------
55 // Implementation of the TPC Kr cluster class
61 // Origin: Adam Matyja, INP PAN, adam.matyja@ifj.edu.pl
64 //-----------------------------------------------------------------
70 #include "AliTPCclusterKr.h"
73 #include "AliCluster.h"
76 #include "AliTPCvtpr.h"
79 #include "TObjArray.h"
94 ClassImp(AliTPCclusterKr)
103 AliTPCclusterKr::AliTPCclusterKr()
164 // default constructor
170 fCluster=new TObjArray();
179 AliTPCclusterKr::AliTPCclusterKr(const AliTPCclusterKr ¶m)
246 fADCcluster = param.fADCcluster;
252 fNPads = param.fNPads;
255 fNRows = param.fNRows;
261 // fCluster = param.fCluster;
264 fCenterX = param.fCenterX;
267 fCenterY = param.fCenterY;
270 fCenterT = param.fCenterT;
273 fCluster=new TObjArray(*(param.fCluster));
279 fTimebins1D = param.fTimebins1D;
282 fPads1D = param.fPads1D;
285 fPadRMS = param.fPadRMS;
288 fRowRMS = param.fRowRMS;
291 fTimebinRMS = param.fTimebinRMS;
294 fTimeStamp = param.fTimeStamp;
304 AliTPCclusterKr &AliTPCclusterKr::operator = (const AliTPCclusterKr & param)
313 // assignment operator
317 if (this == ¶m) return (*this);
319 (AliCluster&)(*this) = (AliCluster&)param;
322 fADCcluster = param.fADCcluster;
328 fNPads = param.fNPads;
331 fNRows = param.fNRows;
337 // fCluster=param.fCluster;
340 fCenterX = param.fCenterX;
343 fCenterY = param.fCenterY;
346 fCenterT = param.fCenterT;
352 fCluster=new TObjArray(*(param.fCluster));
358 fTimebins1D = param.fTimebins1D;
361 fPads1D = param.fPads1D;
364 fPadRMS = param.fPadRMS;
367 fRowRMS = param.fRowRMS;
370 fTimebinRMS = param.fTimebinRMS;
373 fTimeStamp = param.fTimeStamp;
386 AliTPCclusterKr::~AliTPCclusterKr()
404 fCluster->SetOwner(kTRUE);
425 ////____________________________________________________________________________
428 void AliTPCclusterKr::SetCenter(){
434 // calculate geometrical center of the cluster
458 for(Int_t iter = 0; iter < fCluster->GetEntriesFast(); ++iter) {
461 AliTPCvtpr *iclus=(AliTPCvtpr *)fCluster->At(iter);
467 //for( std::vector<AliTPCvtpr*>::iterator iclus = fCluster.begin();
470 //iclus != fCluster.end(); ++iclus ) {
473 adc = (iclus)->GetAdc();
479 rX += ((iclus)->GetX() * adc);
482 rY += ((iclus)->GetY() * adc);
485 rT += ((iclus)->GetT() * adc);
491 fCenterX=rX/fADCcluster;
494 fCenterY=rY/fADCcluster;
497 fCenterT=rT/fADCcluster;
513 void AliTPCclusterKr::SetPadRMS(){
519 // calculate RMS in pad direction
525 // TH1F *histo= new TH1F("","",200,0,200);
528 TArrayI *array= new TArrayI(fCluster->GetEntriesFast());
531 for(Int_t i=0;i<fCluster->GetEntriesFast();i++)
537 array->SetAt(((AliTPCvtpr *)(fCluster->At(i)))->GetPad(),i);
540 //histo->Fill( ((AliTPCvtpr *)(fCluster->At(i)))->GetPad() );
546 // fPadRMS=histo->GetRMS();
549 fPadRMS=TMath::RMS(array->GetSize(),array->GetArray());
567 void AliTPCclusterKr::SetRowRMS(){
573 // calculate RMS in row direction
579 TArrayI *array= new TArrayI(fCluster->GetEntriesFast());
582 // TH1F *histo= new TH1F("","",120,0,120);
585 for(Int_t i=0;i<fCluster->GetEntriesFast();i++)
591 array->SetAt(((AliTPCvtpr *)(fCluster->At(i)))->GetRow(),i);
594 // histo->Fill( ((AliTPCvtpr *)(fCluster->At(i)))->GetRow() );
600 // fRowRMS=histo->GetRMS();
603 fRowRMS=TMath::RMS(array->GetSize(),array->GetArray());
621 void AliTPCclusterKr::SetTimebinRMS(){
627 // calculate RMS in timebin direction
633 TArrayI *array= new TArrayI(fCluster->GetEntriesFast());
636 // TH1F *histo= new TH1F("","",1000,0,1000);
639 for(Int_t i=0;i<fCluster->GetEntriesFast();i++)
645 array->SetAt(((AliTPCvtpr *)(fCluster->At(i)))->GetTime(),i);
648 // histo->Fill( ((AliTPCvtpr *)(fCluster->At(i)))->GetTime() );
654 fTimebinRMS=TMath::RMS(array->GetSize(),array->GetArray());
675 void AliTPCclusterKr::SetRMS(){
681 // calculate RMS in pad,row,timebin direction
687 TArrayI *arrayPad = new TArrayI(fCluster->GetEntriesFast());
690 TArrayI *arrayRow = new TArrayI(fCluster->GetEntriesFast());
693 TArrayI *arrayTime= new TArrayI(fCluster->GetEntriesFast());
696 // TH1F *histoPad= new TH1F("p","p",200,0,200);
699 // TH1F *histoRow= new TH1F("r","r",120,0,120);
702 // TH1F *histoTime= new TH1F("t","t",1000,0,1000);
705 for(Int_t i=0;i<fCluster->GetEntriesFast();i++)
711 arrayPad->SetAt(((AliTPCvtpr *)(fCluster->At(i)))->GetPad(),i);
714 arrayRow->SetAt(((AliTPCvtpr *)(fCluster->At(i)))->GetRow(),i);
717 arrayTime->SetAt(((AliTPCvtpr *)(fCluster->At(i)))->GetTime(),i);
723 //histoPad->Fill( ((AliTPCvtpr *)(fCluster->At(i)))->GetPad() );
726 //histoRow->Fill( ((AliTPCvtpr *)(fCluster->At(i)))->GetRow() );
729 //histoTime->Fill( ((AliTPCvtpr *)(fCluster->At(i)))->GetTime() );
735 // fPadRMS=histoPad->GetRMS();
738 fPadRMS=TMath::RMS(arrayPad->GetSize(),arrayPad->GetArray());
741 fRowRMS=TMath::RMS(arrayRow->GetSize(),arrayRow->GetArray());
744 //histoRow->GetRMS();
747 fTimebinRMS=TMath::RMS(arrayTime->GetSize(),arrayTime->GetArray());
750 //histoTime->GetRMS();
789 void AliTPCclusterKr::Set1D(){
804 Short_t minTime=1000;
816 for(Int_t i=0;i<fCluster->GetEntriesFast();i++)
822 if(((AliTPCvtpr *)(fCluster->At(i)))->GetPad()>maxPad)maxPad =((AliTPCvtpr *)(fCluster->At(i)))->GetPad();
825 if(((AliTPCvtpr *)(fCluster->At(i)))->GetPad()<minPad)minPad =((AliTPCvtpr *)(fCluster->At(i)))->GetPad();
828 if(((AliTPCvtpr *)(fCluster->At(i)))->GetTime()>maxTime)maxTime=((AliTPCvtpr *)(fCluster->At(i)))->GetTime();
831 if(((AliTPCvtpr *)(fCluster->At(i)))->GetTime()<minTime)minTime=((AliTPCvtpr *)(fCluster->At(i)))->GetTime();
837 fPads1D=maxPad-minPad+1;
840 fTimebins1D=maxTime-minTime+1;