]> git.uio.no Git - u/mrichter/AliRoot.git/commitdiff
Adding cut on isolations
authormarian <marian@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 24 Oct 2008 16:52:23 +0000 (16:52 +0000)
committermarian <marian@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 24 Oct 2008 16:52:23 +0000 (16:52 +0000)
Adding the cluster shape
(Adam Matyja)

TPC/AliTPCclusterKr.cxx
TPC/AliTPCclusterKr.h
TPC/AliTPCclustererKr.cxx
TPC/AliTPCclustererKr.h
TPC/FindKrClusters.C
TPC/FindKrClustersRaw.C

index 921c538f78d4fd023902a7d9839cf9fa740840db..33e0e292b920564589aecae03c6400dd6749ac17 100644 (file)
-/**************************************************************************
- * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
- *                                                                        *
- * Author: The ALICE Off-line Project.                                    *
- * Contributors are mentioned in the code where appropriate.              *
- *                                                                        *
- * Permission to use, copy, modify and distribute this software and its   *
- * documentation strictly for non-commercial purposes is hereby granted   *
- * without fee, provided that the above copyright notice appears in all   *
- * copies and that both the copyright notice and this permission notice   *
- * appear in the supporting documentation. The authors make no claims     *
- * about the suitability of this software for any purpose. It is          *
- * provided "as is" without express or implied warranty.                  *
- **************************************************************************/
-
-/* $Id: AliTPCclusterKr.cxx,v 1.7 2008/01/22 17:24:53 matyja Exp $ */
-
-//-----------------------------------------------------------------
-//           Implementation of the TPC Kr cluster class
-//
-// Origin: Adam Matyja, INP PAN, adam.matyja@ifj.edu.pl
-//-----------------------------------------------------------------
-
-#include "AliTPCclusterKr.h"
-#include "AliCluster.h"
-#include "AliTPCvtpr.h"
-#include "TObjArray.h"
-
-ClassImp(AliTPCclusterKr)
-
-
-AliTPCclusterKr::AliTPCclusterKr()
-:AliCluster(),
- fMax(),
- fADCcluster(0),
- fSec(0),
- fNPads(0),
- fNRows(0),
- fSize(0),
- fCenterX(0),
- fCenterY(0),
- fCenterT(0),
- fCluster(0)
-{
-//
-// default constructor
-//
-  fCluster=new TObjArray();
-}
-
-AliTPCclusterKr::AliTPCclusterKr(const AliTPCclusterKr &param)
-:AliCluster(param),
- fMax(),
- fADCcluster(0),
- fSec(0),
- fNPads(0),
- fNRows(0),
- fSize(0),
- fCenterX(0),
- fCenterY(0),
- fCenterT(0),
- fCluster(0)
-{
-//
-// copy constructor
-//
-  fADCcluster = param.fADCcluster;
-  fSec  = param.fSec ;
-  fNPads = param.fNPads;
-  fNRows = param.fNRows;
-  fMax = param.fMax;
-  //  fCluster = param.fCluster;
-  fCenterX = param.fCenterX;
-  fCenterY = param.fCenterY;
-  fCenterT = param.fCenterT;
-  fCluster=new TObjArray(*(param.fCluster));
-  fSize = param.fSize;
-} 
-
-AliTPCclusterKr &AliTPCclusterKr::operator = (const AliTPCclusterKr & param)
-{
-  //
-  // assignment operator
-  // 
-  (AliCluster&)(*this) = (AliCluster&)param;
-  fADCcluster = param.fADCcluster;
-  fSec  = param.fSec ;
-  fNPads = param.fNPads;
-  fNRows = param.fNRows;
-  fMax = param.fMax;
-  //  fCluster=param.fCluster;
-  fCenterX = param.fCenterX;
-  fCenterY = param.fCenterY;
-  fCenterT = param.fCenterT;
-  delete fCluster;
-  fCluster=new TObjArray(*(param.fCluster));
-  fSize=param.fSize;
-  return (*this);
-}
-
-AliTPCclusterKr::~AliTPCclusterKr()
-{
-  //
-  // destructor
-  //
-  if(fCluster) {
-    fCluster->SetOwner(kTRUE);
-    fCluster->Delete();
-    delete fCluster;
-  }
-  fCluster=0;
-}
-
-////____________________________________________________________________________
-void AliTPCclusterKr::SetCenter(){
-  //
-  // calculate geometrical center of the cluster
-  //
-  Double_t rX=0;
-  Double_t rY=0;
-  Double_t rT=0;
-
-  Short_t adc;
-  fADCcluster=0;
-  for(Int_t iter = 0; iter < fCluster->GetEntriesFast(); ++iter) {
-    AliTPCvtpr *iclus=(AliTPCvtpr *)fCluster->At(iter);
-
-    //for( std::vector<AliTPCvtpr*>::iterator iclus  = fCluster.begin();
-    //iclus != fCluster.end(); ++iclus ) {
-    adc = (iclus)->GetAdc();
-    fADCcluster+=adc;
-    rX += ((iclus)->GetX() * adc);
-    rY += ((iclus)->GetY() * adc);
-    rT += ((iclus)->GetT() * adc);
-  }
-  fCenterX=rX/fADCcluster;
-  fCenterY=rY/fADCcluster;
-  fCenterT=rT/fADCcluster;
-
-  return;
-}
+/**************************************************************************\r
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *\r
+ *                                                                        *\r
+ * Author: The ALICE Off-line Project.                                    *\r
+ * Contributors are mentioned in the code where appropriate.              *\r
+ *                                                                        *\r
+ * Permission to use, copy, modify and distribute this software and its   *\r
+ * documentation strictly for non-commercial purposes is hereby granted   *\r
+ * without fee, provided that the above copyright notice appears in all   *\r
+ * copies and that both the copyright notice and this permission notice   *\r
+ * appear in the supporting documentation. The authors make no claims     *\r
+ * about the suitability of this software for any purpose. It is          *\r
+ * provided "as is" without express or implied warranty.                  *\r
+ **************************************************************************/\r
+\r
+/* $Id: AliTPCclusterKr.cxx,v 1.7 2008/01/22 17:24:53 matyja Exp $ */\r
+\r
+//-----------------------------------------------------------------\r
+//           Implementation of the TPC Kr cluster class\r
+//\r
+// Origin: Adam Matyja, INP PAN, adam.matyja@ifj.edu.pl\r
+//-----------------------------------------------------------------\r
+\r
+#include "AliTPCclusterKr.h"\r
+#include "AliCluster.h"\r
+#include "AliTPCvtpr.h"\r
+#include "TObjArray.h"\r
+//#include "TH1F.h"\r
+#include "TMath.h"\r
+#include "TArrayI.h"\r
+\r
+ClassImp(AliTPCclusterKr)\r
+\r
+\r
+AliTPCclusterKr::AliTPCclusterKr()\r
+:AliCluster(),\r
+ fMax(),\r
+ fADCcluster(0),\r
+ fSec(0),\r
+ fNPads(0),\r
+ fNRows(0),\r
+ fTimebins1D(0),\r
+ fPads1D(0),\r
+ fPadRMS(0),\r
+ fRowRMS(0),\r
+ fTimebinRMS(0),\r
+ fSize(0),\r
+ fCenterX(0),\r
+ fCenterY(0),\r
+ fCenterT(0),\r
+ fCluster(0)\r
+{\r
+//\r
+// default constructor\r
+//\r
+  fCluster=new TObjArray();\r
+}\r
+\r
+AliTPCclusterKr::AliTPCclusterKr(const AliTPCclusterKr &param)\r
+:AliCluster(param),\r
+ fMax(),\r
+ fADCcluster(0),\r
+ fSec(0),\r
+ fNPads(0),\r
+ fNRows(0),\r
+ fTimebins1D(0),\r
+ fPads1D(0),\r
+ fPadRMS(0),\r
+ fRowRMS(0),\r
+ fTimebinRMS(0),\r
+ fSize(0),\r
+ fCenterX(0),\r
+ fCenterY(0),\r
+ fCenterT(0),\r
+ fCluster(0)\r
+{\r
+//\r
+// copy constructor\r
+//\r
+  fADCcluster = param.fADCcluster;\r
+  fSec  = param.fSec ;\r
+  fNPads = param.fNPads;\r
+  fNRows = param.fNRows;\r
+  fMax = param.fMax;\r
+  //  fCluster = param.fCluster;\r
+  fCenterX = param.fCenterX;\r
+  fCenterY = param.fCenterY;\r
+  fCenterT = param.fCenterT;\r
+  fCluster=new TObjArray(*(param.fCluster));\r
+  fSize = param.fSize;\r
+  fTimebins1D = param.fTimebins1D;\r
+  fPads1D = param.fPads1D;\r
+  fPadRMS = param.fPadRMS;\r
+  fRowRMS = param.fRowRMS;\r
+  fTimebinRMS = param.fTimebinRMS;\r
+} \r
+\r
+AliTPCclusterKr &AliTPCclusterKr::operator = (const AliTPCclusterKr & param)\r
+{\r
+  //\r
+  // assignment operator\r
+  // \r
+  (AliCluster&)(*this) = (AliCluster&)param;\r
+  fADCcluster = param.fADCcluster;\r
+  fSec  = param.fSec ;\r
+  fNPads = param.fNPads;\r
+  fNRows = param.fNRows;\r
+  fMax = param.fMax;\r
+  //  fCluster=param.fCluster;\r
+  fCenterX = param.fCenterX;\r
+  fCenterY = param.fCenterY;\r
+  fCenterT = param.fCenterT;\r
+  delete fCluster;\r
+  fCluster=new TObjArray(*(param.fCluster));\r
+  fSize=param.fSize;\r
+  fTimebins1D = param.fTimebins1D;\r
+  fPads1D = param.fPads1D;\r
+  fPadRMS = param.fPadRMS;\r
+  fRowRMS = param.fRowRMS;\r
+  fTimebinRMS = param.fTimebinRMS;\r
+  return (*this);\r
+}\r
+\r
+AliTPCclusterKr::~AliTPCclusterKr()\r
+{\r
+  //\r
+  // destructor\r
+  //\r
+  if(fCluster) {\r
+    fCluster->SetOwner(kTRUE);\r
+    fCluster->Delete();\r
+    delete fCluster;\r
+  }\r
+  fCluster=0;\r
+}\r
+\r
+////____________________________________________________________________________\r
+void AliTPCclusterKr::SetCenter(){\r
+  //\r
+  // calculate geometrical center of the cluster\r
+  //\r
+  Double_t rX=0;\r
+  Double_t rY=0;\r
+  Double_t rT=0;\r
+\r
+  Short_t adc;\r
+  fADCcluster=0;\r
+  for(Int_t iter = 0; iter < fCluster->GetEntriesFast(); ++iter) {\r
+    AliTPCvtpr *iclus=(AliTPCvtpr *)fCluster->At(iter);\r
+\r
+    //for( std::vector<AliTPCvtpr*>::iterator iclus  = fCluster.begin();\r
+    //iclus != fCluster.end(); ++iclus ) {\r
+    adc = (iclus)->GetAdc();\r
+    fADCcluster+=adc;\r
+    rX += ((iclus)->GetX() * adc);\r
+    rY += ((iclus)->GetY() * adc);\r
+    rT += ((iclus)->GetT() * adc);\r
+  }\r
+  fCenterX=rX/fADCcluster;\r
+  fCenterY=rY/fADCcluster;\r
+  fCenterT=rT/fADCcluster;\r
+\r
+  return;\r
+}\r
+\r
+void AliTPCclusterKr::SetPadRMS(){\r
+  //\r
+  // calculate RMS in pad direction\r
+  //\r
+  //  TH1F *histo= new TH1F("","",200,0,200);\r
+  TArrayI *array= new TArrayI(fCluster->GetEntriesFast());\r
+  for(Int_t i=0;i<fCluster->GetEntriesFast();i++)\r
+    {\r
+      array->SetAt(((AliTPCvtpr *)(fCluster->At(i)))->GetPad(),i);\r
+      //histo->Fill( ((AliTPCvtpr *)(fCluster->At(i)))->GetPad() );\r
+    }\r
+  //  fPadRMS=histo->GetRMS();\r
+  fPadRMS=TMath::RMS(array->GetSize(),array->GetArray());\r
+  //  delete histo;\r
+  delete array;\r
+  return;\r
+}\r
+\r
+void AliTPCclusterKr::SetRowRMS(){\r
+  //\r
+  // calculate RMS in row direction\r
+  //\r
+  TArrayI *array= new TArrayI(fCluster->GetEntriesFast());\r
+  //  TH1F *histo= new TH1F("","",120,0,120);\r
+  for(Int_t i=0;i<fCluster->GetEntriesFast();i++)\r
+    {\r
+      array->SetAt(((AliTPCvtpr *)(fCluster->At(i)))->GetRow(),i);\r
+      //      histo->Fill( ((AliTPCvtpr *)(fCluster->At(i)))->GetRow() );\r
+    }\r
+  //  fRowRMS=histo->GetRMS();\r
+  fRowRMS=TMath::RMS(array->GetSize(),array->GetArray());\r
+  //  delete histo;\r
+  delete array;\r
+  return;\r
+}\r
+\r
+void AliTPCclusterKr::SetTimebinRMS(){\r
+  //\r
+  // calculate RMS in timebin direction\r
+  //\r
+  TArrayI *array= new TArrayI(fCluster->GetEntriesFast());\r
+  //  TH1F *histo= new TH1F("","",1000,0,1000);\r
+  for(Int_t i=0;i<fCluster->GetEntriesFast();i++)\r
+    {\r
+      array->SetAt(((AliTPCvtpr *)(fCluster->At(i)))->GetTime(),i);\r
+      //      histo->Fill( ((AliTPCvtpr *)(fCluster->At(i)))->GetTime() );\r
+    }\r
+  fTimebinRMS=TMath::RMS(array->GetSize(),array->GetArray());\r
+  //histo->GetRMS();\r
+  //  delete histo;\r
+  delete array;\r
+  return;\r
+}\r
+\r
+void AliTPCclusterKr::SetRMS(){\r
+  //\r
+  // calculate RMS in pad,row,timebin direction\r
+  //\r
+  TArrayI *arrayPad = new TArrayI(fCluster->GetEntriesFast());\r
+  TArrayI *arrayRow = new TArrayI(fCluster->GetEntriesFast());\r
+  TArrayI *arrayTime= new TArrayI(fCluster->GetEntriesFast());\r
+  //  TH1F *histoPad= new TH1F("p","p",200,0,200);\r
+  //  TH1F *histoRow= new TH1F("r","r",120,0,120);\r
+  //  TH1F *histoTime= new TH1F("t","t",1000,0,1000);\r
+  for(Int_t i=0;i<fCluster->GetEntriesFast();i++)\r
+    {\r
+      arrayPad->SetAt(((AliTPCvtpr *)(fCluster->At(i)))->GetPad(),i);\r
+      arrayRow->SetAt(((AliTPCvtpr *)(fCluster->At(i)))->GetRow(),i);\r
+      arrayTime->SetAt(((AliTPCvtpr *)(fCluster->At(i)))->GetTime(),i);\r
+\r
+      //histoPad->Fill( ((AliTPCvtpr *)(fCluster->At(i)))->GetPad() );\r
+      //histoRow->Fill( ((AliTPCvtpr *)(fCluster->At(i)))->GetRow() );\r
+      //histoTime->Fill( ((AliTPCvtpr *)(fCluster->At(i)))->GetTime() );\r
+    }\r
+  //  fPadRMS=histoPad->GetRMS();\r
+  fPadRMS=TMath::RMS(arrayPad->GetSize(),arrayPad->GetArray());\r
+  fRowRMS=TMath::RMS(arrayRow->GetSize(),arrayRow->GetArray());\r
+    //histoRow->GetRMS();\r
+  fTimebinRMS=TMath::RMS(arrayTime->GetSize(),arrayTime->GetArray());\r
+    //histoTime->GetRMS();\r
+\r
+  delete arrayPad;\r
+  delete arrayRow;\r
+  delete arrayTime;\r
+  //  delete histoPad;\r
+  //  delete histoRow;\r
+  //  delete histoTime;\r
+\r
+  return;\r
+}\r
+\r
+\r
+void AliTPCclusterKr::Set1D(){\r
+  //\r
+  //\r
+  //\r
+  Short_t maxTime=0;\r
+  Short_t minTime=1000;\r
+  Short_t maxPad=0;\r
+  Short_t minPad=1000;\r
\r
+  for(Int_t i=0;i<fCluster->GetEntriesFast();i++)\r
+    {\r
+      if(((AliTPCvtpr *)(fCluster->At(i)))->GetPad()>maxPad)maxPad   =((AliTPCvtpr *)(fCluster->At(i)))->GetPad();\r
+      if(((AliTPCvtpr *)(fCluster->At(i)))->GetPad()<minPad)minPad   =((AliTPCvtpr *)(fCluster->At(i)))->GetPad();\r
+      if(((AliTPCvtpr *)(fCluster->At(i)))->GetTime()>maxTime)maxTime=((AliTPCvtpr *)(fCluster->At(i)))->GetTime();\r
+      if(((AliTPCvtpr *)(fCluster->At(i)))->GetTime()<minTime)minTime=((AliTPCvtpr *)(fCluster->At(i)))->GetTime();\r
+    }\r
+  fPads1D=maxPad-minPad+1;\r
+  fTimebins1D=maxTime-minTime+1;\r
+  return;\r
+}\r
index c5945f80493e880f96a32b6d0930577c42823814..281ba71f031ec7734ce7e18d02401247c56eda0d 100644 (file)
-#ifndef ALITPCCLUSTERKR_H
-#define ALITPCCLUSTERKR_H
-/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
- * See cxx source for full Copyright notice                               */
-
-/* $Id: AliTPCclusterKr.h,v 1.8 2008/01/22 16:07:15 matyja Exp $ */
-
-//-------------------------------------------------------
-//                    TPC Kr Cluster Class
-//
-//   Origin: Adam Matyja, INP PAN, adam.matyja@ifj.edu.pl
-//-------------------------------------------------------
-
-#include "AliCluster.h"
-#include "TObjArray.h"
-#include "AliTPCvtpr.h"
-
-
-//_____________________________________________________________________________
-class AliTPCclusterKr: public AliCluster{
-public:
-  AliTPCclusterKr();
-  AliTPCclusterKr(const AliTPCclusterKr & param);//copy constructor
-  AliTPCclusterKr &operator = (const AliTPCclusterKr & param); 
-  virtual ~AliTPCclusterKr();
-
-  virtual void SetCenter();//set center of the cluster weighted by charge
-
-  virtual void SetMax(AliTPCvtpr q){fMax=q;}//set values of max. in cluster
-  virtual void SetADCcluster(Int_t q){fADCcluster=q;}
-  virtual void SetSec(Short_t q){fSec=q;}
-  virtual void SetNPads(Short_t q){fNPads=q;}
-  virtual void SetNRows(Short_t q){fNRows=q;}
-  virtual void SetSize(){fSize=fCluster->GetEntriesFast();}
-  virtual void SetCenterX(Double_t q){fCenterX=q;}
-  virtual void SetCenterY(Double_t q){fCenterY=q;}
-  virtual void SetCenterT(Double_t q){fCenterT=q;}
-  //void AddDigitToCluster(AliTPCvtpr *q){fCluster.push_back(q);}
-  virtual void AddDigitToCluster(AliTPCvtpr *q){
-    fCluster->AddLast(q);
-    //fCluster->Compress();
-  }
-
-  AliTPCvtpr GetMax() const {return fMax;}
-  Int_t GetADCcluster() const {return  fADCcluster;}
-  Short_t GetSec() const {return fSec;}
-  Short_t GetNPads() const {return fNPads;}
-  Short_t GetNRows() const {return fNRows;}
-  Short_t GetSize() const {return fSize;}
-  Double_t GetCenterX() const {return fCenterX;}
-  Double_t GetCenterY() const {return fCenterY;}
-  Double_t GetCenterT() const {return fCenterT;}
-  AliTPCvtpr *GetDigitFromCluster(Int_t i) const {return (AliTPCvtpr*)fCluster->At(i);}
-
-private:
-  AliTPCvtpr fMax;//max (ADC,timebin,pad,row) in cluster
-  Int_t fADCcluster; //ADC of cluster
-  Short_t fSec;  //sector of the cluster
-  Short_t fNPads; //number of pads in cluster
-  Short_t fNRows; //number of rows in cluster
-  Short_t fSize; //size of vector
-  Double_t fCenterX;// X coordinate of the cluster center in cm
-  Double_t fCenterY;// Y coordinate of the cluster center in cm
-  Double_t fCenterT;// time coordinate of the cluster center in timebins
-  //std::vector< AliTPCvtpr*> fCluster;//cluster contents(adc,nt,np,nr)
-  TObjArray *fCluster;//cluster contents(adc,nt,np,nr)
-
-
-  ClassDef(AliTPCclusterKr,4)  // Time Projection Chamber Kr clusters
-};
-
-
-#endif
-
-
+#ifndef ALITPCCLUSTERKR_H\r
+#define ALITPCCLUSTERKR_H\r
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *\r
+ * See cxx source for full Copyright notice                               */\r
+\r
+/* $Id: AliTPCclusterKr.h,v 1.8 2008/01/22 16:07:15 matyja Exp $ */\r
+\r
+//-------------------------------------------------------\r
+//                    TPC Kr Cluster Class\r
+//\r
+//   Origin: Adam Matyja, INP PAN, adam.matyja@ifj.edu.pl\r
+//-------------------------------------------------------\r
+\r
+#include "AliCluster.h"\r
+#include "TObjArray.h"\r
+#include "AliTPCvtpr.h"\r
+//#include "TH1F.h"\r
+#include "TMath.h"\r
+#include "TArrayI.h"\r
+\r
+//_____________________________________________________________________________\r
+class AliTPCclusterKr: public AliCluster{\r
+public:\r
+  AliTPCclusterKr();\r
+  AliTPCclusterKr(const AliTPCclusterKr & param);//copy constructor\r
+  AliTPCclusterKr &operator = (const AliTPCclusterKr & param); \r
+  virtual ~AliTPCclusterKr();\r
+\r
+  virtual void SetCenter();//set center of the cluster weighted by charge\r
+\r
+  virtual void SetMax(AliTPCvtpr q){fMax=q;}//set values of max. in cluster\r
+  virtual void SetADCcluster(Int_t q){fADCcluster=q;}\r
+  virtual void SetSec(Short_t q){fSec=q;}\r
+  virtual void SetNPads(Short_t q){fNPads=q;}\r
+  virtual void SetNRows(Short_t q){fNRows=q;}\r
+  virtual void SetSize(){fSize=fCluster->GetEntriesFast();}\r
+  virtual void SetCenterX(Double_t q){fCenterX=q;}\r
+  virtual void SetCenterY(Double_t q){fCenterY=q;}\r
+  virtual void SetCenterT(Double_t q){fCenterT=q;}\r
+\r
+  virtual void SetTimebins1D(Short_t q){fTimebins1D=q;}\r
+  virtual void SetPads1D(Short_t q){fPads1D=q;}\r
+  virtual void Set1D();\r
+  virtual void SetPadRMS(Double_t q){fPadRMS=q;}\r
+  virtual void SetRowRMS(Double_t q){fRowRMS=q;}\r
+  virtual void SetTimebinRMS(Double_t q){fTimebinRMS=q;}\r
+  virtual void SetPadRMS();\r
+  virtual void SetRowRMS();\r
+  virtual void SetTimebinRMS();\r
+  virtual void SetRMS();\r
+  //void AddDigitToCluster(AliTPCvtpr *q){fCluster.push_back(q);}\r
+  virtual void AddDigitToCluster(AliTPCvtpr *q){\r
+    fCluster->AddLast(q);\r
+    //fCluster->Compress();\r
+  }\r
+\r
+  AliTPCvtpr GetMax() const {return fMax;}\r
+  Int_t GetADCcluster() const {return  fADCcluster;}\r
+  Short_t GetSec() const {return fSec;}\r
+  Short_t GetNPads() const {return fNPads;}\r
+  Short_t GetNRows() const {return fNRows;}\r
+  Short_t GetSize() const {return fSize;}\r
+\r
+  Short_t GetTimebins1D(){return fTimebins1D;}\r
+  Short_t GetPads1D(){return fPads1D;}\r
+  Double_t GetPadRMS(){return fPadRMS;}\r
+  Double_t GetRowRMS(){return fRowRMS;}\r
+  Double_t GetTimebinRMS(){return fTimebinRMS;}\r
+\r
+  Double_t GetCenterX() const {return fCenterX;}\r
+  Double_t GetCenterY() const {return fCenterY;}\r
+  Double_t GetCenterT() const {return fCenterT;}\r
+  AliTPCvtpr *GetDigitFromCluster(Int_t i) const {return (AliTPCvtpr*)fCluster->At(i);}\r
+\r
+private:\r
+  AliTPCvtpr fMax;//max (ADC,timebin,pad,row) in cluster\r
+  Int_t fADCcluster; //ADC of cluster\r
+  Short_t fSec;  //sector of the cluster\r
+  Short_t fNPads; //number of pads in cluster\r
+  Short_t fNRows; //number of rows in cluster or row max - min\r
+\r
+  Short_t fTimebins1D; //Timebin max - min\r
+  Short_t fPads1D; //Pad max - min\r
+  Double_t fPadRMS; //Pad RMS\r
+  Double_t fRowRMS; //Row RMS\r
+  Double_t fTimebinRMS; //Timebin RMS\r
+\r
+  Short_t fSize; //size of vector\r
+  Double_t fCenterX;// X coordinate of the cluster center in cm\r
+  Double_t fCenterY;// Y coordinate of the cluster center in cm\r
+  Double_t fCenterT;// time coordinate of the cluster center in timebins\r
+  //std::vector< AliTPCvtpr*> fCluster;//cluster contents(adc,nt,np,nr)\r
+  TObjArray *fCluster;//cluster contents(adc,nt,np,nr)\r
+\r
+\r
+  ClassDef(AliTPCclusterKr,5)  // Time Projection Chamber Kr clusters\r
+};\r
+\r
+\r
+#endif\r
+\r
+\r
index bac6dd9a0ea3d5c309de978d490a182d13858343..b093b6a79a08aaa213535ac9eecdc3b471a3eb2d 100644 (file)
-/**************************************************************************
- * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
- *                                                                        *
- * Author: The ALICE Off-line Project.                                    *
- * Contributors are mentioned in the code where appropriate.              *
- *                                                                        *
- * Permission to use, copy, modify and distribute this software and its   *
- * documentation strictly for non-commercial purposes is hereby granted   *
- * without fee, provided that the above copyright notice appears in all   *
- * copies and that both the copyright notice and this permission notice   *
- * appear in the supporting documentation. The authors make no claims     *
- * about the suitability of this software for any purpose. It is          *
- * provided "as is" without express or implied warranty.                  *
- **************************************************************************/
-
-/* $Id: AliTPCclustererKr.cxx,v 1.7 2008/02/06 17:24:53 matyja Exp $ */
-
-//-----------------------------------------------------------------
-//           Implementation of the TPC Kr cluster class
-//
-// Origin: Adam Matyja, INP PAN, adam.matyja@ifj.edu.pl
-//-----------------------------------------------------------------
-
-/*
-Instruction - how to use that:
-There are two macros prepared. One is for preparing clusters from MC 
-samples:  FindKrClusters.C. The output is kept in TPC.RecPoints.root.
-The other macro is prepared for data analysis: FindKrClustersRaw.C. 
-The output is created for each processed file in root file named adc.root. 
-For each data subsample the same named file is created. So be careful 
-do not overwrite them. 
-
-Additional selection criteria to select the GOLD cluster
-Example:
-// open  file with clusters
-TFile f("Krypton.root");
-TTree * tree = (TTree*)f.Get("Kr")
-TCut cutR0("cutR0","fADCcluster/fSize<100");        // adjust it according v seetings - 
-TCut cutR1("cutR1","fADCcluster/fSize>7");          // cosmic tracks and noise removal
-TCut cutR2("cutR2","fMax.fAdc/fADCcluster<0.2");    // digital noise removal
-TCut cutR3("cutR3","fMax.fAdc/fADCcluster>0.01");   // noise removal
-TCut cutS1("cutS1","fSize<200");    // adjust it according v seetings - cosmic tracks
-TCut cutAll = cutR0+cutR1+cutR2+cutR3+cutS1;
-This values are typical values to be applied in selectors
-
-
-*
-**** MC ****
-*
-
-To run clusterizaton for MC type:
-.x FindKrClusters.C
-
-If you don't want to use the standard selection criteria then you 
-have to do following:
-
-// load the standard setup
-AliRunLoader* rl = AliRunLoader::Open("galice.root");
-AliTPCLoader *tpcl = (AliTPCLoader*)rl->GetLoader("TPCLoader");
-tpcl->LoadDigits();
-rl->LoadgAlice();
-gAlice=rl->GetAliRun();
-TDirectory *cwd = gDirectory;
-AliTPCv4 *tpc = (AliTPCv4*)gAlice->GetDetector("TPC");
-Int_t ver = tpc->IsVersion();
-rl->CdGAFile();
-AliTPCParam *param=(AliTPCParamSR *)gDirectory->Get("75x40_100x60_150x60");
-AliTPCDigitsArray *digarr=new AliTPCDigitsArray;
-digarr->Setup(param);
-cwd->cd();
-
-//loop over events
-Int_t nevmax=rl->GetNumberOfEvents();
-for(Int_t nev=0;nev<nevmax ;nev++){
-  rl->GetEvent(nev);
-  TTree* input_tree= tpcl->TreeD();//load tree with digits
-  digarr->ConnectTree(input_tree);
-  TTree *output_tree =tpcl->TreeR();//load output tree
-
-  AliTPCclustererKr *clusters = new AliTPCclustererKr();
-  clusters->SetParam(param);
-  clusters->SetInput(input_tree);
-  clusters->SetOutput(output_tree);
-  clusters->SetDigArr(digarr);
-  
-//If you want to change the cluster finder parameters for MC there are 
-//several of them:
-
-//1. signal threshold (everything below the given number is treated as 0)
-  clusters->SetMinAdc(3);
-
-//2. number of neighbouring timebins to be considered
-  clusters->SetMinTimeBins(2);
-
-//3. distance of the cluster center to the center of a pad in pad-padrow plane 
-//(in cm). Remenber that this is still quantified by pad size.
-  clusters->SetMaxPadRangeCm(2.5);
-
-//4. distance of the cluster center to the center of a padrow in pad-padrow 
-//plane (in cm). Remenber that this is still quantified by pad size.
-  clusters->SetMaxRowRangeCm(3.5);
-
-//5. distance of the cluster center to the max time bin on a pad (in tackts)
-//ie. fabs(centerT - time)<7
-  clusters->SetMaxTimeRange(7);
-
-//6. cut reduce peak at 0. There are noises which appear mostly as two 
-//timebins on one pad.
-  clusters->SetValueToSize(3.5);
-
-
-  clusters->finderIO();
-  tpcl->WriteRecPoints("OVERWRITE");
-}
-delete rl;//cleans everything
-
-*
-********* DATA *********
-*
-
-To run clusterizaton for DATA for file named raw_data.root type:
-.x FindKrClustersRaw.C("raw_data.root")
-
-If you want to change some criteria do the following:
-
-//
-// remove Altro warnings
-//
-AliLog::SetClassDebugLevel("AliTPCRawStream",-5);
-AliLog::SetClassDebugLevel("AliRawReaderDate",-5);
-AliLog::SetClassDebugLevel("AliTPCAltroMapping",-5);
-AliLog::SetModuleDebugLevel("RAW",-5);
-
-//
-// Get database with noises
-//
-//  char *ocdbpath = gSystem->Getenv("OCDB_PATH");
-char *ocdbpath ="local:///afs/cern.ch/alice/tpctest/OCDB";
-if (ocdbpath==0){
-ocdbpath="alien://folder=/alice/data/2007/LHC07w/OCDB/";
-}
-AliCDBManager * man = AliCDBManager::Instance();
-man->SetDefaultStorage(ocdbpath);
-man->SetRun(0);
-AliTPCCalPad * noiseTPC = AliTPCcalibDB::Instance()->GetPadNoise();
-AliTPCAltroMapping** mapping =AliTPCcalibDB::Instance()->GetMapping();
-
-//define tree
-TFile *hfile=new TFile("adc.root","RECREATE","ADC file");
-// Create a ROOT Tree
-TTree *mytree = new TTree("Kr","Krypton cluster tree");
-
-//define infput file
-const char *fileName="data.root";
-AliRawReader *reader = new AliRawReaderRoot(fileName);
-//AliRawReader *reader = new AliRawReaderDate(fileName);
-reader->Reset();
-AliAltroRawStreamFast* stream = new AliAltroRawStreamFast(reader);
-stream->SelectRawData("TPC");
-
-//one general output
-AliTPCclustererKr *clusters = new AliTPCclustererKr();
-clusters->SetOutput(mytree);
-clusters->SetRecoParam(0);//standard reco parameters
-AliTPCParamSR *param=new AliTPCParamSR();
-clusters->SetParam(param);//TPC parameters(sectors, timebins, etc.)
-
-//set cluster finder parameters (from data):
-//1. zero suppression parameter
-  clusters->SetZeroSup(param->GetZeroSup());
-
-//2. first bin
-  clusters->SetFirstBin(60);
-
-//3. last bin
-  clusters->SetLastBin(950);
-
-//4. maximal noise
-  clusters->SetMaxNoiseAbs(2);
-
-//5. maximal amount of sigma of noise
-  clusters->SetMaxNoiseSigma(3);
-
-//The remaining parameters are the same paramters as for MC (see MC section 
-//points 1-6)
-  clusters->SetMinAdc(3);
-  clusters->SetMinTimeBins(2);
-  clusters->SetMaxPadRangeCm(2.5);
-  clusters->SetMaxRowRangeCm(3.5);
-  clusters->SetMaxTimeRange(7);
-  clusters->SetValueToSize(3.5);
-
-while (reader->NextEvent()) {
-  clusters->FinderIO(reader);
-}
-
-hfile->Write();
-hfile->Close();
-delete stream;
-
-
-*/
-
-#include "AliTPCclustererKr.h"
-#include "AliTPCclusterKr.h"
-//#include <vector>
-#include <list>
-#include "TObject.h"
-#include "AliPadMax.h"
-#include "AliSimDigits.h"
-#include "AliTPCv4.h"
-#include "AliTPCParam.h"
-#include "AliTPCDigitsArray.h"
-#include "AliTPCvtpr.h"
-#include "AliTPCClustersRow.h"
-#include "TTree.h"
-#include "TH1F.h"
-#include "TH2F.h"
-#include "TTreeStream.h"
-
-#include "AliTPCTransform.h"
-
-//used in raw data finder
-#include "AliTPCROC.h"
-#include "AliTPCCalPad.h"
-#include "AliTPCAltroMapping.h"
-#include "AliTPCcalibDB.h"
-#include "AliTPCRawStream.h"
-#include "AliTPCRecoParam.h"
-#include "AliTPCReconstructor.h"
-#include "AliRawReader.h"
-#include "AliTPCCalROC.h"
-
-ClassImp(AliTPCclustererKr)
-
-
-AliTPCclustererKr::AliTPCclustererKr()
-  :TObject(),
-  fRawData(kFALSE),
-  fInput(0),
-  fOutput(0),
-  fParam(0),
-  fDigarr(0),
-  fRecoParam(0),
-  fZeroSup(2),
-  fFirstBin(60),
-  fLastBin(950),
-  fMaxNoiseAbs(2),
-  fMaxNoiseSigma(3),
-  fMinAdc(3),
-  fMinTimeBins(2),
-//  fMaxPadRange(4),
-//  fMaxRowRange(3),
-  fMaxTimeRange(7),
-  fValueToSize(3.5),
-  fMaxPadRangeCm(2.5),
-  fMaxRowRangeCm(3.5),
-  fDebugLevel(-1),
-  fHistoRow(0),
-  fHistoPad(0),
-  fHistoTime(0),
-  fHistoRowPad(0)
-{
-//
-// default constructor
-//
-}
-
-AliTPCclustererKr::AliTPCclustererKr(const AliTPCclustererKr &param)
-  :TObject(),
-  fRawData(kFALSE),
-  fInput(0),
-  fOutput(0),
-  fParam(0),
-  fDigarr(0),
-  fRecoParam(0),
-  fZeroSup(2),
-  fFirstBin(60),
-  fLastBin(950),
-  fMaxNoiseAbs(2),
-  fMaxNoiseSigma(3),
-  fMinAdc(3),
-  fMinTimeBins(2),
-//  fMaxPadRange(4),
-//  fMaxRowRange(3),
-  fMaxTimeRange(7),
-  fValueToSize(3.5),
-  fMaxPadRangeCm(2.5),
-  fMaxRowRangeCm(3.5),
-  fDebugLevel(-1),
-  fHistoRow(0),
-  fHistoPad(0),
-  fHistoTime(0),
-  fHistoRowPad(0)
-{
-//
-// copy constructor
-//
-  fParam = param.fParam;
-  fRecoParam = param.fRecoParam;
-  fRawData = param.fRawData;
-  fInput  = param.fInput ;
-  fOutput = param.fOutput;
-  fDigarr = param.fDigarr;
-  fZeroSup       = param.fZeroSup       ;
-  fFirstBin     = param.fFirstBin      ;
-  fLastBin      = param.fLastBin       ;
-  fMaxNoiseAbs  = param.fMaxNoiseAbs   ;
-  fMaxNoiseSigma = param.fMaxNoiseSigma ;
-  fMinAdc = param.fMinAdc;
-  fMinTimeBins = param.fMinTimeBins;
-//  fMaxPadRange  = param.fMaxPadRange ;
-//  fMaxRowRange  = param.fMaxRowRange ;
-  fMaxTimeRange = param.fMaxTimeRange;
-  fValueToSize  = param.fValueToSize;
-  fMaxPadRangeCm = param.fMaxPadRangeCm;
-  fMaxRowRangeCm = param.fMaxRowRangeCm;
-  fDebugLevel = param.fDebugLevel;
-  fHistoRow    = param.fHistoRow   ;
-  fHistoPad    = param.fHistoPad  ;
-  fHistoTime   = param.fHistoTime;
-  fHistoRowPad = param.fHistoRowPad;
-
-} 
-
-AliTPCclustererKr & AliTPCclustererKr::operator = (const AliTPCclustererKr & param)
-{
-  //
-  // assignment operator
-  //
-  fParam = param.fParam;
-  fRecoParam = param.fRecoParam;
-  fRawData = param.fRawData;
-  fInput  = param.fInput ;
-  fOutput = param.fOutput;
-  fDigarr = param.fDigarr;
-  fZeroSup       = param.fZeroSup       ;
-  fFirstBin     = param.fFirstBin      ;
-  fLastBin      = param.fLastBin       ;
-  fMaxNoiseAbs  = param.fMaxNoiseAbs   ;
-  fMaxNoiseSigma = param.fMaxNoiseSigma ;
-  fMinAdc = param.fMinAdc;
-  fMinTimeBins = param.fMinTimeBins;
-//  fMaxPadRange  = param.fMaxPadRange ;
-//  fMaxRowRange  = param.fMaxRowRange ;
-  fMaxTimeRange = param.fMaxTimeRange;
-  fValueToSize  = param.fValueToSize;
-  fMaxPadRangeCm = param.fMaxPadRangeCm;
-  fMaxRowRangeCm = param.fMaxRowRangeCm;
-  fDebugLevel = param.fDebugLevel;
-  fHistoRow    = param.fHistoRow   ;
-  fHistoPad    = param.fHistoPad  ;
-  fHistoTime   = param.fHistoTime;
-  fHistoRowPad = param.fHistoRowPad;
-  return (*this);
-}
-
-AliTPCclustererKr::~AliTPCclustererKr()
-{
-  //
-  // destructor
-  //
-  delete fOutput;
-}
-
-void AliTPCclustererKr::SetRecoParam(AliTPCRecoParam *recoParam)
-{
-  //
-  // set reconstruction parameters
-  //
-  if (recoParam) {
-    fRecoParam = recoParam;
-  }else{
-    //set default parameters if not specified
-    fRecoParam = AliTPCReconstructor::GetRecoParam();
-    if (!fRecoParam)  fRecoParam = AliTPCRecoParam::GetLowFluxParam();
-  }
-  return;
-}
-
-
-////____________________________________________________________________________
-////       I/O
-void AliTPCclustererKr::SetInput(TTree * tree)
-{
-  //
-  // set input tree with digits
-  //
-  fInput = tree;  
-  if  (!fInput->GetBranch("Segment")){
-    cerr<<"AliTPCclusterKr::FindClusterKr(): no proper input tree !\n";
-    fInput=0;
-    return;
-  }
-}
-
-void AliTPCclustererKr::SetOutput(TTree * /*tree*/) 
-{
-  //
-  //dummy
-  //
-  fOutput = new TTreeSRedirector("Krypton.root");
-}
-
-////____________________________________________________________________________
-//// with new I/O
-Int_t AliTPCclustererKr::FinderIO()
-{
-  // Krypton cluster finder for simulated events from MC
-
-  if (!fInput) { 
-    Error("Digits2Clusters", "input tree not initialised");
-    return 10;
-  }
-  
-  if (!fOutput) {
-    Error("Digits2Clusters", "output tree not initialised");
-    return 11;
-  }
-
-  FindClusterKrIO();
-  return 0;
-}
-
-Int_t AliTPCclustererKr::FinderIO(AliRawReader* rawReader)
-{
-  // Krypton cluster finder for the TPC raw data
-  //
-  // fParam must be defined before
-
-  if(rawReader)fRawData=kTRUE; //set flag to data
-
-  if (!fOutput) {
-    Error("Digits2Clusters", "output tree not initialised");
-    return 11;
-  }
-
-  fParam->SetMaxTBin(fRecoParam->GetLastBin());//set number of timebins from reco -> param
-  //   used later for memory allocation
-
-  Bool_t isAltro=kFALSE;
-
-  AliTPCROC * roc = AliTPCROC::Instance();
-  AliTPCCalPad * noiseTPC = AliTPCcalibDB::Instance()->GetPadNoise();
-  AliTPCAltroMapping** mapping =AliTPCcalibDB::Instance()->GetMapping();
-  //
-  AliTPCRawStream input(rawReader,(AliAltroMapping**)mapping);
-
-  const Int_t kNIS = fParam->GetNInnerSector();//number of inner sectors
-  const Int_t kNOS = fParam->GetNOuterSector();//number of outer sectors
-  const Int_t kNS = kNIS + kNOS;//all sectors
-
-
-  //crate TPC view
-  AliTPCDigitsArray *digarr=new AliTPCDigitsArray(kFALSE);//data not sim
-  digarr->Setup(fParam);//as usually parameters
-
-  //
-  // Loop over sectors
-  //
-  for(Int_t iSec = 0; iSec < kNS; iSec++) {
-    AliTPCCalROC * noiseROC;
-    if(noiseTPC==0x0){
-      noiseROC =new AliTPCCalROC(iSec);//noise=0
-    }
-    else{
-      noiseROC = noiseTPC->GetCalROC(iSec);  // noise per given sector
-    }
-    Int_t nRows = 0; //number of rows in sector
-    Int_t nDDLs = 0; //number of DDLs
-    Int_t indexDDL = 0; //DDL index
-    if (iSec < kNIS) {
-      nRows = fParam->GetNRowLow();
-      nDDLs = 2;
-      indexDDL = iSec * 2;
-    }else {
-      nRows = fParam->GetNRowUp();
-      nDDLs = 4;
-      indexDDL = (iSec-kNIS) * 4 + kNIS * 2;
-    }
-
-    //
-    // Load the raw data for corresponding DDLs
-    //
-    rawReader->Reset();
-    rawReader->Select("TPC",indexDDL,indexDDL+nDDLs-1);
-
-    if(input.Next()) {
-      isAltro=kTRUE;
-      // Allocate memory for rows in sector (pads(depends on row) x timebins)
-      for(Int_t iRow = 0; iRow < nRows; iRow++) {
-       digarr->CreateRow(iSec,iRow); 
-      }//end loop over rows
-    }
-    rawReader->Reset();
-    rawReader->Select("TPC",indexDDL,indexDDL+nDDLs-1);
-
-    //
-    // Begin loop over altro data
-    //
-    while (input.Next()) {
-      
-      //check sector consistency
-      if (input.GetSector() != iSec)
-       AliFatal(Form("Sector index mismatch ! Expected (%d), but got (%d) !",iSec,input.GetSector()));
-      
-      Int_t iRow = input.GetRow();
-      Int_t iPad = input.GetPad();
-      Int_t iTimeBin = input.GetTime();
-
-      //
-      if(fDebugLevel==72){
-       fHistoRow->Fill(iRow);
-       fHistoPad->Fill(iPad);
-       fHistoTime->Fill(iTimeBin);
-       fHistoRowPad->Fill(iPad,iRow);
-      }else if(fDebugLevel>=0&&fDebugLevel<72){
-       if(iSec==fDebugLevel){
-         fHistoRow->Fill(iRow);
-         fHistoPad->Fill(iPad);
-         fHistoTime->Fill(iTimeBin);
-         fHistoRowPad->Fill(iPad,iRow);
-       }
-      }else if(fDebugLevel==73){
-       if(iSec<36){
-         fHistoRow->Fill(iRow);
-         fHistoPad->Fill(iPad);
-         fHistoTime->Fill(iTimeBin);
-         fHistoRowPad->Fill(iPad,iRow);
-       }
-      }else if(fDebugLevel==74){
-       if(iSec>=36){
-         fHistoRow->Fill(iRow);
-         fHistoPad->Fill(iPad);
-         fHistoTime->Fill(iTimeBin);
-         fHistoRowPad->Fill(iPad,iRow);
-       }
-      }
-
-      //check row consistency
-      if (iRow < 0 ) continue;
-      if (iRow < 0 || iRow >= nRows){
-       AliError(Form("Pad-row index (%d) outside the range (%d -> %d) !",
-                     iRow, 0, nRows -1));
-       continue;
-      }
-
-      //check pad consistency
-      if (iPad < 0 || iPad >= (Int_t)(roc->GetNPads(iSec,iRow))) {
-       AliError(Form("Pad index (%d) outside the range (%d -> %d) !",
-                     iPad, 0, roc->GetNPads(iSec,iRow) ));
-       continue;
-      }
-
-      //check time consistency
-      if ( iTimeBin < fRecoParam->GetFirstBin() || iTimeBin >= fRecoParam->GetLastBin()){
-       //cout<<iTimeBin<<endl;
-       continue;
-       AliFatal(Form("Timebin index (%d) outside the range (%d -> %d) !",
-                     iTimeBin, 0, fRecoParam->GetLastBin() -1));
-      }
-
-      //signal
-      Int_t signal = input.GetSignal();
-      if (signal <= fZeroSup || 
-         iTimeBin < fFirstBin ||
-         iTimeBin > fLastBin
-         ) {
-       digarr->GetRow(iSec,iRow)->SetDigitFast(0,iTimeBin,iPad);
-       continue;
-      }
-      if (!noiseROC) continue;
-      Double_t noiseOnPad = noiseROC->GetValue(iRow,iPad);//noise on given pad and row in sector
-      if (noiseOnPad > fMaxNoiseAbs){
-       digarr->GetRow(iSec,iRow)->SetDigitFast(0,iTimeBin,iPad);
-       continue; // consider noisy pad as dead
-      }
-      if(signal <= fMaxNoiseSigma * noiseOnPad){
-       digarr->GetRow(iSec,iRow)->SetDigitFast(0,iTimeBin,iPad);
-       continue;
-      }
-      digarr->GetRow(iSec,iRow)->SetDigitFast(signal,iTimeBin,iPad);      
-    }//end of loop over altro data
-  }//end of loop over sectors
-  
-  SetDigArr(digarr);
-  if(isAltro) FindClusterKrIO();
-  delete digarr;
-
-  return 0;
-}
-
-
-Int_t AliTPCclustererKr::CleanSector(Int_t sector){
-  //
-  // clean isolated digits
-  //  
-  const Int_t kNRows=fParam->GetNRow(sector);//number of rows in sector
-  for(Int_t iRow=0; iRow<kNRows; ++iRow){
-    AliSimDigits *digrow;
-    if(fRawData){
-      digrow = (AliSimDigits*)fDigarr->GetRow(sector,iRow);//real data
-    }else{
-      digrow = (AliSimDigits*)fDigarr->LoadRow(sector,iRow);//MC
-    }
-    if(!digrow) continue;
-    digrow->ExpandBuffer(); //decrunch
-    const Int_t kNPads = digrow->GetNCols();  // number of pads
-    const Int_t kNTime = digrow->GetNRows(); // number of timebins
-    for(Int_t iPad=1;iPad<kNPads-1;iPad++){
-      Short_t*  val = digrow->GetDigitsColumn(iPad);
-
-      for(Int_t iTimeBin=1;iTimeBin<kNTime-1;iTimeBin++){
-       if (val[iTimeBin]<=0) continue;
-       if (val[iTimeBin-1]+val[iTimeBin+1]<fZeroSup) {val[iTimeBin]=0; continue;}
-       if (val[iTimeBin-kNTime]+val[iTimeBin+kNTime]<fZeroSup) {val[iTimeBin]=0; continue;}
-       //
-       if (val[iTimeBin-1-kNTime]+val[iTimeBin+1+kNTime]<fZeroSup) {val[iTimeBin]=0; continue;}
-       if (val[iTimeBin+1-kNTime]+val[iTimeBin-1+kNTime]<fZeroSup) {val[iTimeBin]=0; continue;}                
-      }
-    }
-  }
-  return 0;
-}
-////____________________________________________________________________________
-Int_t AliTPCclustererKr::FindClusterKrIO()
-{
-
-  //
-  //fParam and  fDigarr must be set to run this method
-  //
-
-  Int_t clusterCounter=0;
-  const Int_t nTotalSector=fParam->GetNSector();//number of sectors
-  for(Int_t iSec=0; iSec<nTotalSector; ++iSec){
-    CleanSector(iSec);
-    //vector of maxima for each sector
-    //std::vector<AliPadMax*> maximaInSector;
-    TObjArray *maximaInSector=new TObjArray();//to store AliPadMax*
-
-    //
-    //  looking for the maxima on the pad
-    //
-
-
-
-    const Int_t kNRows=fParam->GetNRow(iSec);//number of rows in sector
-    for(Int_t iRow=0; iRow<kNRows; ++iRow){
-      AliSimDigits *digrow;
-      if(fRawData){
-       digrow = (AliSimDigits*)fDigarr->GetRow(iSec,iRow);//real data
-      }else{
-       digrow = (AliSimDigits*)fDigarr->LoadRow(iSec,iRow);//MC
-      }
-      if(digrow){//if pointer exist
-       digrow->ExpandBuffer(); //decrunch
-       const Int_t kNPads = digrow->GetNCols();  // number of pads
-       const Int_t kNTime = digrow->GetNRows(); // number of timebins
-       for(Int_t iPad=0;iPad<kNPads;iPad++){
-         
-         Int_t timeBinMax=-1;//timebin of maximum 
-         Int_t valueMaximum=-1;//value of maximum in adc
-         Int_t increaseBegin=-1;//timebin when increase starts
-         Int_t sumAdc=0;//sum of adc on the pad in maximum surrounding
-         bool ifIncreaseBegin=true;//flag - check if increasing started
-         bool ifMaximum=false;//flag - check if it could be maximum
-         Short_t* val = digrow->GetDigitsColumn(iPad);
-         for(Int_t iTimeBin=1;iTimeBin<kNTime-1;iTimeBin++){
-           if (!ifMaximum)  {
-             if (val[iTimeBin]==-1) break;   // 0 until the end
-             for( ; iTimeBin<kNTime-2&&val[iTimeBin]<fMinAdc ;iTimeBin++){}
-           }
-           //
-           Short_t adc = val[iTimeBin];
-           if(adc<fMinAdc){//standard was 3
-             if(ifMaximum){
-               if(iTimeBin-increaseBegin<fMinTimeBins){//at least 2 time bins
-                 timeBinMax=-1;
-                 valueMaximum=-1;
-                 increaseBegin=-1;
-                 sumAdc=0;
-                 ifIncreaseBegin=true;
-                 ifMaximum=false;
-                 continue;
-               }
-               //insert maximum, default values and set flags
-               //Double_t xCord,yCord;
-               //GetXY(iSec,iRow,iPad,xCord,yCord);
-               Double_t x[]={iRow,iPad,iTimeBin};
-               Int_t i[]={iSec};
-               AliTPCTransform trafo;
-               trafo.Transform(x,i,0,1);
-               
-               AliPadMax *oneMaximum = new AliPadMax(AliTPCvtpr(valueMaximum,
-                                                                timeBinMax,
-                                                                iPad,
-                                                                iRow,
-                                                                x[0],//xCord,
-                                                                x[1],//yCord,
-                                                                x[2]/*timeBinMax*/),
-                                                     increaseBegin,
-                                                     iTimeBin-1,
-                                                     sumAdc);
-               maximaInSector->AddLast(oneMaximum);
-               
-               timeBinMax=-1;
-               valueMaximum=-1;
-               increaseBegin=-1;
-               sumAdc=0;
-               ifIncreaseBegin=true;
-               ifMaximum=false;
-             }
-             continue;
-           }
-           
-           if(ifIncreaseBegin){
-             ifIncreaseBegin=false;
-             increaseBegin=iTimeBin;
-           }
-           
-           if(adc>valueMaximum){
-             timeBinMax=iTimeBin;
-             valueMaximum=adc;
-             ifMaximum=true;
-           }
-           sumAdc+=adc;
-           if(iTimeBin==kNTime-1 && ifMaximum && kNTime-increaseBegin>fMinTimeBins){//on the edge
-             //at least 3 timebins
-             //insert maximum, default values and set flags
-             //Double_t xCord,yCord;
-             //GetXY(iSec,iRow,iPad,xCord,yCord);
-             Double_t x[]={iRow,iPad,iTimeBin};
-             Int_t i[]={iSec};
-             AliTPCTransform trafo;
-             trafo.Transform(x,i,0,1);
-             AliPadMax *oneMaximum = new AliPadMax(AliTPCvtpr(valueMaximum,
-                                                              timeBinMax,
-                                                              iPad,
-                                                              iRow,
-                                                              x[0],//xCord,
-                                                              x[1],//yCord,
-                                                              x[2]/*timeBinMax*/),
-                                                   increaseBegin,
-                                                   iTimeBin-1,
-                                                   sumAdc);
-             maximaInSector->AddLast(oneMaximum);
-               
-             timeBinMax=-1;
-             valueMaximum=-1;
-             increaseBegin=-1;
-             sumAdc=0;
-             ifIncreaseBegin=true;
-             ifMaximum=false;
-             continue;
-           }
-           
-         }//end loop over timebins
-       }//end loop over pads
-//      }else{
-//     cout<<"Pointer does not exist!!"<<endl;
-      }//end if poiner exists
-    }//end loop over rows
-
-    MakeClusters(maximaInSector,iSec,clusterCounter);
-    //
-    maximaInSector->SetOwner(kTRUE);
-    maximaInSector->Delete();
-    delete maximaInSector;
-  }//end sector for
-  cout<<"Number of clusters in event: "<<clusterCounter<<endl;
-  return 0;
-}
-
-void AliTPCclustererKr::MakeClusters(TObjArray * maximaInSector, Int_t iSec, Int_t &clusterCounter){
-  //
-  // Make clusters
-  //
-
-  Int_t maxDig=0;
-  Int_t maxSumAdc=0;
-  Int_t maxTimeBin=0;
-  Int_t maxPad=0;
-  Int_t maxRow=0;
-  Double_t maxX=0;
-  Double_t maxY=0;
-  Double_t maxT=0;
-  Int_t entriesArr = maximaInSector->GetEntriesFast();
-  for(Int_t it1 = 0; it1 < entriesArr; ++it1 ) {
-    
-    AliPadMax *mp1=(AliPadMax *)maximaInSector->UncheckedAt(it1);
-    if (!mp1) continue;
-    AliTPCclusterKr clusterKr;
-    
-    Int_t nUsedPads=1;
-    Int_t clusterValue=0;
-    clusterValue+=(mp1)->GetSum();
-    list<Int_t> nUsedRows;
-    nUsedRows.push_back((mp1)->GetRow());
-    
-    maxDig      =(mp1)->GetAdc() ;
-    maxSumAdc   =(mp1)->GetSum() ;
-    maxTimeBin  =(mp1)->GetTime();
-    maxPad      =(mp1)->GetPad() ;
-    maxRow      =(mp1)->GetRow() ;
-    maxX        =(mp1)->GetX();
-    maxY        =(mp1)->GetY();
-    maxT        =(mp1)->GetT();
-    
-    AliSimDigits *digrowTmp;
-    if(fRawData){
-      digrowTmp = (AliSimDigits*)fDigarr->GetRow(iSec,(mp1)->GetRow());
-    }else{
-      digrowTmp = (AliSimDigits*)fDigarr->LoadRow(iSec,(mp1)->GetRow());
-    }
-    
-    digrowTmp->ExpandBuffer(); //decrunch
-    
-    for(Int_t itb=(mp1)->GetBegin(); itb<((mp1)->GetEnd())+1; itb++){
-      Int_t adcTmp = digrowTmp->GetDigitUnchecked(itb,(mp1)->GetPad());
-      AliTPCvtpr *vtpr=new AliTPCvtpr(adcTmp,itb,(mp1)->GetPad(),(mp1)->GetRow(),(mp1)->GetX(),(mp1)->GetY(),(mp1)->GetT());
-      clusterKr.AddDigitToCluster(vtpr);
-    }
-    clusterKr.SetCenter();//set centr of the cluster
-    
-    for(Int_t it2 = it1+1; it2 < entriesArr; ++it2 ) {
-      AliPadMax *mp2=(AliPadMax *)maximaInSector->UncheckedAt(it2);
-      if (!mp2) continue;
-      if (TMath::Abs(clusterKr.GetCenterX() - (mp2)->GetX())>fMaxPadRangeCm) continue;
-      if (TMath::Abs(clusterKr.GetCenterY() - (mp2)->GetY()) > fMaxRowRangeCm) continue;      
-      if (TMath::Abs(clusterKr.GetCenterT() - (mp2)->GetT()) >fMaxTimeRange) continue;
-
-      {
-       clusterValue+=(mp2)->GetSum();
-       
-       nUsedPads++;
-       nUsedRows.push_back((mp2)->GetRow());
-       
-       AliSimDigits *digrowTmp1;
-       if(fRawData){
-         digrowTmp1 = (AliSimDigits*)fDigarr->GetRow(iSec,(mp2)->GetRow());
-       }else{
-         digrowTmp1 = (AliSimDigits*)fDigarr->LoadRow(iSec,(mp2)->GetRow());
-       }
-       digrowTmp1->ExpandBuffer(); //decrunch
-       
-       for(Int_t itb=(mp2)->GetBegin(); itb<(mp2)->GetEnd()+1; itb++){
-         Int_t adcTmp = digrowTmp1->GetDigitUnchecked(itb,(mp2)->GetPad());
-         AliTPCvtpr *vtpr=new AliTPCvtpr(adcTmp,itb,(mp2)->GetPad(),(mp2)->GetRow(),(mp2)->GetX(),(mp2)->GetY(),(mp2)->GetT());
-         clusterKr.AddDigitToCluster(vtpr);
-       }
-       
-       clusterKr.SetCenter();//set center of the cluster
-       
-       //which one is bigger
-       if( (mp2)->GetAdc() > maxDig ){
-         maxDig      =(mp2)->GetAdc() ;
-         maxSumAdc   =(mp2)->GetSum() ;
-         maxTimeBin  =(mp2)->GetTime();
-         maxPad      =(mp2)->GetPad() ;
-         maxRow      =(mp2)->GetRow() ;
-         maxX        =(mp2)->GetX() ;
-         maxY        =(mp2)->GetY() ;
-         maxT        =(mp2)->GetT() ;
-       } else if ( (mp2)->GetAdc() == maxDig ){
-         if( (mp2)->GetSum() > maxSumAdc){
-           maxDig      =(mp2)->GetAdc() ;
-           maxSumAdc   =(mp2)->GetSum() ;
-           maxTimeBin  =(mp2)->GetTime();
-           maxPad      =(mp2)->GetPad() ;
-           maxRow      =(mp2)->GetRow() ;
-           maxX        =(mp2)->GetX() ;
-           maxY        =(mp2)->GetY() ;
-           maxT        =(mp2)->GetT() ;
-         }
-       }
-       delete maximaInSector->RemoveAt(it2);
-      }
-    }//inside loop
-    delete maximaInSector->RemoveAt(it1);          
-    clusterKr.SetSize();
-    //through out clusters on the edge and noise
-    //if(clusterValue/clusterKr.fCluster.size()<fValueToSize)continue;
-    if(clusterValue/(clusterKr.GetSize())<fValueToSize)continue;
-    
-    clusterKr.SetADCcluster(clusterValue);
-    clusterKr.SetNPads(nUsedPads);
-    clusterKr.SetMax(AliTPCvtpr(maxDig,maxTimeBin,maxPad,maxRow,maxX,maxY,maxT));
-    clusterKr.SetSec(iSec);
-    clusterKr.SetSize();
-    
-    nUsedRows.sort();
-    nUsedRows.unique();
-    clusterKr.SetNRows(nUsedRows.size());
-    clusterKr.SetCenter();
-    
-    clusterCounter++;
-    
-    
-    //save each cluster into file
-    if (fOutput){
-      (*fOutput)<<"Kr"<<
-       "Cl.="<<&clusterKr<<
-       "\n";
-    }
-    //end of save each cluster into file adc.root
-  }//outer loop
-}
-
-
-
-////____________________________________________________________________________
-
-
-void AliTPCclustererKr::GetXY(Int_t sec,Int_t row,Int_t pad,Double_t& xGlob,Double_t& yGlob){
-  //
-  //gives global XY coordinate of the pad
-  //
-
-  Double_t yLocal = fParam->GetPadRowRadii(sec,row);//radius of row in sector in cm
-
-  Int_t padmax = fParam->GetNPads(sec,row);//number of pads in a given row
-  Float_t padXSize;
-  if(sec<fParam->GetNInnerSector())padXSize=0.4;
-  else padXSize=0.6;
-  Double_t xLocal=(pad+0.5-padmax/2.)*padXSize;//x-value of the center of pad
-
-  Float_t sin,cos;
-  fParam->AdjustCosSin((Int_t)sec,cos,sin);//return sinus and cosinus of the sector
-
-  Double_t xGlob1 =  xLocal * cos + yLocal * sin;
-  Double_t yGlob1 = -xLocal * sin + yLocal * cos;
-
-
-  Double_t rot=0;
-  rot=TMath::Pi()/2.;
-
-  xGlob =  xGlob1 * TMath::Cos(rot) + yGlob1 * TMath::Sin(rot);
-  yGlob = -xGlob1 * TMath::Sin(rot) + yGlob1 * TMath::Cos(rot);
-
-   yGlob=-1*yGlob;
-   if(sec<18||(sec>=36&&sec<54)) xGlob =-1*xGlob;
-
-
-  return;
-}
+/**************************************************************************\r
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *\r
+ *                                                                        *\r
+ * Author: The ALICE Off-line Project.                                    *\r
+ * Contributors are mentioned in the code where appropriate.              *\r
+ *                                                                        *\r
+ * Permission to use, copy, modify and distribute this software and its   *\r
+ * documentation strictly for non-commercial purposes is hereby granted   *\r
+ * without fee, provided that the above copyright notice appears in all   *\r
+ * copies and that both the copyright notice and this permission notice   *\r
+ * appear in the supporting documentation. The authors make no claims     *\r
+ * about the suitability of this software for any purpose. It is          *\r
+ * provided "as is" without express or implied warranty.                  *\r
+ **************************************************************************/\r
+\r
+/* $Id: AliTPCclustererKr.cxx,v 1.7 2008/02/06 17:24:53 matyja Exp $ */\r
+\r
+//-----------------------------------------------------------------\r
+//           Implementation of the TPC Kr cluster class\r
+//\r
+// Origin: Adam Matyja, INP PAN, adam.matyja@ifj.edu.pl\r
+//-----------------------------------------------------------------\r
+\r
+/*\r
+Instruction - how to use that:\r
+There are two macros prepared. One is for preparing clusters from MC \r
+samples:  FindKrClusters.C. The output is kept in TPC.RecPoints.root.\r
+The other macro is prepared for data analysis: FindKrClustersRaw.C. \r
+The output is created for each processed file in root file named adc.root. \r
+For each data subsample the same named file is created. So be careful \r
+do not overwrite them. \r
+\r
+Additional selection criteria to select the GOLD cluster\r
+Example:\r
+// open  file with clusters\r
+TFile f("Krypton.root");\r
+TTree * tree = (TTree*)f.Get("Kr")\r
+TCut cutR0("cutR0","fADCcluster/fSize<100");        // adjust it according v seetings - \r
+TCut cutR1("cutR1","fADCcluster/fSize>7");          // cosmic tracks and noise removal\r
+TCut cutR2("cutR2","fMax.fAdc/fADCcluster<0.2");    // digital noise removal\r
+TCut cutR3("cutR3","fMax.fAdc/fADCcluster>0.01");   // noise removal\r
+TCut cutS1("cutS1","fSize<200");    // adjust it according v seetings - cosmic tracks\r
+TCut cutAll = cutR0+cutR1+cutR2+cutR3+cutS1;\r
+This values are typical values to be applied in selectors\r
+\r
+\r
+*\r
+**** MC ****\r
+*\r
+\r
+To run clusterizaton for MC type:\r
+.x FindKrClusters.C\r
+\r
+If you don't want to use the standard selection criteria then you \r
+have to do following:\r
+\r
+// load the standard setup\r
+AliRunLoader* rl = AliRunLoader::Open("galice.root");\r
+AliTPCLoader *tpcl = (AliTPCLoader*)rl->GetLoader("TPCLoader");\r
+tpcl->LoadDigits();\r
+rl->LoadgAlice();\r
+gAlice=rl->GetAliRun();\r
+TDirectory *cwd = gDirectory;\r
+AliTPCv4 *tpc = (AliTPCv4*)gAlice->GetDetector("TPC");\r
+Int_t ver = tpc->IsVersion();\r
+rl->CdGAFile();\r
+AliTPCParam *param=(AliTPCParamSR *)gDirectory->Get("75x40_100x60_150x60");\r
+AliTPCDigitsArray *digarr=new AliTPCDigitsArray;\r
+digarr->Setup(param);\r
+cwd->cd();\r
+\r
+//loop over events\r
+Int_t nevmax=rl->GetNumberOfEvents();\r
+for(Int_t nev=0;nev<nevmax ;nev++){\r
+  rl->GetEvent(nev);\r
+  TTree* input_tree= tpcl->TreeD();//load tree with digits\r
+  digarr->ConnectTree(input_tree);\r
+  TTree *output_tree =tpcl->TreeR();//load output tree\r
+\r
+  AliTPCclustererKr *clusters = new AliTPCclustererKr();\r
+  clusters->SetParam(param);\r
+  clusters->SetInput(input_tree);\r
+  clusters->SetOutput(output_tree);\r
+  clusters->SetDigArr(digarr);\r
+  \r
+//If you want to change the cluster finder parameters for MC there are \r
+//several of them:\r
+\r
+//1. signal threshold (everything below the given number is treated as 0)\r
+  clusters->SetMinAdc(3);\r
+\r
+//2. number of neighbouring timebins to be considered\r
+  clusters->SetMinTimeBins(2);\r
+\r
+//3. distance of the cluster center to the center of a pad in pad-padrow plane \r
+//(in cm). Remenber that this is still quantified by pad size.\r
+  clusters->SetMaxPadRangeCm(2.5);\r
+\r
+//4. distance of the cluster center to the center of a padrow in pad-padrow \r
+//plane (in cm). Remenber that this is still quantified by pad size.\r
+  clusters->SetMaxRowRangeCm(3.5);\r
+\r
+//5. distance of the cluster center to the max time bin on a pad (in tackts)\r
+//ie. fabs(centerT - time)<7\r
+  clusters->SetMaxTimeRange(7);\r
+\r
+//6. cut reduce peak at 0. There are noises which appear mostly as two \r
+//timebins on one pad.\r
+  clusters->SetValueToSize(3.5);\r
+\r
+\r
+  clusters->finderIO();\r
+  tpcl->WriteRecPoints("OVERWRITE");\r
+}\r
+delete rl;//cleans everything\r
+\r
+*\r
+********* DATA *********\r
+*\r
+\r
+To run clusterizaton for DATA for file named raw_data.root type:\r
+.x FindKrClustersRaw.C("raw_data.root")\r
+\r
+If you want to change some criteria do the following:\r
+\r
+//\r
+// remove Altro warnings\r
+//\r
+AliLog::SetClassDebugLevel("AliTPCRawStream",-5);\r
+AliLog::SetClassDebugLevel("AliRawReaderDate",-5);\r
+AliLog::SetClassDebugLevel("AliTPCAltroMapping",-5);\r
+AliLog::SetModuleDebugLevel("RAW",-5);\r
+\r
+//\r
+// Get database with noises\r
+//\r
+//  char *ocdbpath = gSystem->Getenv("OCDB_PATH");\r
+char *ocdbpath ="local:///afs/cern.ch/alice/tpctest/OCDB";\r
+if (ocdbpath==0){\r
+ocdbpath="alien://folder=/alice/data/2007/LHC07w/OCDB/";\r
+}\r
+AliCDBManager * man = AliCDBManager::Instance();\r
+man->SetDefaultStorage(ocdbpath);\r
+man->SetRun(0);\r
+AliTPCCalPad * noiseTPC = AliTPCcalibDB::Instance()->GetPadNoise();\r
+AliTPCAltroMapping** mapping =AliTPCcalibDB::Instance()->GetMapping();\r
+\r
+//define tree\r
+TFile *hfile=new TFile("adc.root","RECREATE","ADC file");\r
+// Create a ROOT Tree\r
+TTree *mytree = new TTree("Kr","Krypton cluster tree");\r
+\r
+//define infput file\r
+const char *fileName="data.root";\r
+AliRawReader *reader = new AliRawReaderRoot(fileName);\r
+//AliRawReader *reader = new AliRawReaderDate(fileName);\r
+reader->Reset();\r
+AliAltroRawStreamFast* stream = new AliAltroRawStreamFast(reader);\r
+stream->SelectRawData("TPC");\r
+\r
+//one general output\r
+AliTPCclustererKr *clusters = new AliTPCclustererKr();\r
+clusters->SetOutput(mytree);\r
+clusters->SetRecoParam(0);//standard reco parameters\r
+AliTPCParamSR *param=new AliTPCParamSR();\r
+clusters->SetParam(param);//TPC parameters(sectors, timebins, etc.)\r
+\r
+//set cluster finder parameters (from data):\r
+//1. zero suppression parameter\r
+  clusters->SetZeroSup(param->GetZeroSup());\r
+\r
+//2. first bin\r
+  clusters->SetFirstBin(60);\r
+\r
+//3. last bin\r
+  clusters->SetLastBin(950);\r
+\r
+//4. maximal noise\r
+  clusters->SetMaxNoiseAbs(2);\r
+\r
+//5. maximal amount of sigma of noise\r
+  clusters->SetMaxNoiseSigma(3);\r
+\r
+//The remaining parameters are the same paramters as for MC (see MC section \r
+//points 1-6)\r
+  clusters->SetMinAdc(3);\r
+  clusters->SetMinTimeBins(2);\r
+  clusters->SetMaxPadRangeCm(2.5);\r
+  clusters->SetMaxRowRangeCm(3.5);\r
+  clusters->SetMaxTimeRange(7);\r
+  clusters->SetValueToSize(3.5);\r
+\r
+while (reader->NextEvent()) {\r
+  clusters->FinderIO(reader);\r
+}\r
+\r
+hfile->Write();\r
+hfile->Close();\r
+delete stream;\r
+\r
+\r
+*/\r
+\r
+#include "AliTPCclustererKr.h"\r
+#include "AliTPCclusterKr.h"\r
+//#include <vector>\r
+#include <list>\r
+#include "TObject.h"\r
+#include "AliPadMax.h"\r
+#include "AliSimDigits.h"\r
+#include "AliTPCv4.h"\r
+#include "AliTPCParam.h"\r
+#include "AliTPCDigitsArray.h"\r
+#include "AliTPCvtpr.h"\r
+#include "AliTPCClustersRow.h"\r
+#include "TTree.h"\r
+#include "TH1F.h"\r
+#include "TH2F.h"\r
+#include "TTreeStream.h"\r
+\r
+#include "AliTPCTransform.h"\r
+\r
+//used in raw data finder\r
+#include "AliTPCROC.h"\r
+#include "AliTPCCalPad.h"\r
+#include "AliTPCAltroMapping.h"\r
+#include "AliTPCcalibDB.h"\r
+#include "AliTPCRawStream.h"\r
+#include "AliTPCRecoParam.h"\r
+#include "AliTPCReconstructor.h"\r
+#include "AliRawReader.h"\r
+#include "AliTPCCalROC.h"\r
+\r
+ClassImp(AliTPCclustererKr)\r
+\r
+\r
+AliTPCclustererKr::AliTPCclustererKr()\r
+  :TObject(),\r
+  fRawData(kFALSE),\r
+  fInput(0),\r
+  fOutput(0),\r
+  fParam(0),\r
+  fDigarr(0),\r
+  fRecoParam(0),\r
+  fZeroSup(2),\r
+  fFirstBin(60),\r
+  fLastBin(950),\r
+  fMaxNoiseAbs(2),\r
+  fMaxNoiseSigma(3),\r
+  fMinAdc(3),\r
+  fMinTimeBins(2),\r
+//  fMaxPadRange(4),\r
+//  fMaxRowRange(3),\r
+  fMaxTimeRange(7),\r
+  fValueToSize(3.5),\r
+  fMaxPadRangeCm(2.5),\r
+  fMaxRowRangeCm(3.5),\r
+  fIsolCut(3),\r
+  fDebugLevel(-1),\r
+  fHistoRow(0),\r
+  fHistoPad(0),\r
+  fHistoTime(0),\r
+  fHistoRowPad(0)\r
+{\r
+//\r
+// default constructor\r
+//\r
+}\r
+\r
+AliTPCclustererKr::AliTPCclustererKr(const AliTPCclustererKr &param)\r
+  :TObject(),\r
+  fRawData(kFALSE),\r
+  fInput(0),\r
+  fOutput(0),\r
+  fParam(0),\r
+  fDigarr(0),\r
+  fRecoParam(0),\r
+  fZeroSup(2),\r
+  fFirstBin(60),\r
+  fLastBin(950),\r
+  fMaxNoiseAbs(2),\r
+  fMaxNoiseSigma(3),\r
+  fMinAdc(3),\r
+  fMinTimeBins(2),\r
+//  fMaxPadRange(4),\r
+//  fMaxRowRange(3),\r
+  fMaxTimeRange(7),\r
+  fValueToSize(3.5),\r
+  fMaxPadRangeCm(2.5),\r
+  fMaxRowRangeCm(3.5),\r
+  fIsolCut(3),\r
+  fDebugLevel(-1),\r
+  fHistoRow(0),\r
+  fHistoPad(0),\r
+  fHistoTime(0),\r
+  fHistoRowPad(0)\r
+{\r
+//\r
+// copy constructor\r
+//\r
+  fParam = param.fParam;\r
+  fRecoParam = param.fRecoParam;\r
+  fRawData = param.fRawData;\r
+  fInput  = param.fInput ;\r
+  fOutput = param.fOutput;\r
+  fDigarr = param.fDigarr;\r
+  fZeroSup       = param.fZeroSup       ;\r
+  fFirstBin     = param.fFirstBin      ;\r
+  fLastBin      = param.fLastBin       ;\r
+  fMaxNoiseAbs  = param.fMaxNoiseAbs   ;\r
+  fMaxNoiseSigma = param.fMaxNoiseSigma ;\r
+  fMinAdc = param.fMinAdc;\r
+  fMinTimeBins = param.fMinTimeBins;\r
+//  fMaxPadRange  = param.fMaxPadRange ;\r
+//  fMaxRowRange  = param.fMaxRowRange ;\r
+  fMaxTimeRange = param.fMaxTimeRange;\r
+  fValueToSize  = param.fValueToSize;\r
+  fMaxPadRangeCm = param.fMaxPadRangeCm;\r
+  fMaxRowRangeCm = param.fMaxRowRangeCm;\r
+  fIsolCut = param.fIsolCut;\r
+  fDebugLevel = param.fDebugLevel;\r
+  fHistoRow    = param.fHistoRow   ;\r
+  fHistoPad    = param.fHistoPad  ;\r
+  fHistoTime   = param.fHistoTime;\r
+  fHistoRowPad = param.fHistoRowPad;\r
+\r
+\r
+} \r
+\r
+AliTPCclustererKr & AliTPCclustererKr::operator = (const AliTPCclustererKr & param)\r
+{\r
+  //\r
+  // assignment operator\r
+  //\r
+  fParam = param.fParam;\r
+  fRecoParam = param.fRecoParam;\r
+  fRawData = param.fRawData;\r
+  fInput  = param.fInput ;\r
+  fOutput = param.fOutput;\r
+  fDigarr = param.fDigarr;\r
+  fZeroSup       = param.fZeroSup       ;\r
+  fFirstBin     = param.fFirstBin      ;\r
+  fLastBin      = param.fLastBin       ;\r
+  fMaxNoiseAbs  = param.fMaxNoiseAbs   ;\r
+  fMaxNoiseSigma = param.fMaxNoiseSigma ;\r
+  fMinAdc = param.fMinAdc;\r
+  fMinTimeBins = param.fMinTimeBins;\r
+//  fMaxPadRange  = param.fMaxPadRange ;\r
+//  fMaxRowRange  = param.fMaxRowRange ;\r
+  fMaxTimeRange = param.fMaxTimeRange;\r
+  fValueToSize  = param.fValueToSize;\r
+  fMaxPadRangeCm = param.fMaxPadRangeCm;\r
+  fMaxRowRangeCm = param.fMaxRowRangeCm;\r
+  fIsolCut = param.fIsolCut;\r
+  fDebugLevel = param.fDebugLevel;\r
+  fHistoRow    = param.fHistoRow   ;\r
+  fHistoPad    = param.fHistoPad  ;\r
+  fHistoTime   = param.fHistoTime;\r
+  fHistoRowPad = param.fHistoRowPad;\r
+  return (*this);\r
+}\r
+\r
+AliTPCclustererKr::~AliTPCclustererKr()\r
+{\r
+  //\r
+  // destructor\r
+  //\r
+  delete fOutput;\r
+}\r
+\r
+void AliTPCclustererKr::SetRecoParam(AliTPCRecoParam *recoParam)\r
+{\r
+  //\r
+  // set reconstruction parameters\r
+  //\r
+  if (recoParam) {\r
+    fRecoParam = recoParam;\r
+  }else{\r
+    //set default parameters if not specified\r
+    fRecoParam = AliTPCReconstructor::GetRecoParam();\r
+    if (!fRecoParam)  fRecoParam = AliTPCRecoParam::GetLowFluxParam();\r
+  }\r
+  return;\r
+}\r
+\r
+\r
+////____________________________________________________________________________\r
+////       I/O\r
+void AliTPCclustererKr::SetInput(TTree * tree)\r
+{\r
+  //\r
+  // set input tree with digits\r
+  //\r
+  fInput = tree;  \r
+  if  (!fInput->GetBranch("Segment")){\r
+    cerr<<"AliTPCclusterKr::FindClusterKr(): no proper input tree !\n";\r
+    fInput=0;\r
+    return;\r
+  }\r
+}\r
+\r
+void AliTPCclustererKr::SetOutput(TTree * /*tree*/) \r
+{\r
+  //\r
+  //dummy\r
+  //\r
+  fOutput = new TTreeSRedirector("Krypton.root");\r
+}\r
+\r
+////____________________________________________________________________________\r
+//// with new I/O\r
+Int_t AliTPCclustererKr::FinderIO()\r
+{\r
+  // Krypton cluster finder for simulated events from MC\r
+\r
+  if (!fInput) { \r
+    Error("Digits2Clusters", "input tree not initialised");\r
+    return 10;\r
+  }\r
+  \r
+  if (!fOutput) {\r
+    Error("Digits2Clusters", "output tree not initialised");\r
+    return 11;\r
+  }\r
+\r
+  FindClusterKrIO();\r
+  return 0;\r
+}\r
+\r
+Int_t AliTPCclustererKr::FinderIO(AliRawReader* rawReader)\r
+{\r
+  // Krypton cluster finder for the TPC raw data\r
+  //\r
+  // fParam must be defined before\r
+\r
+  if(rawReader)fRawData=kTRUE; //set flag to data\r
+\r
+  if (!fOutput) {\r
+    Error("Digits2Clusters", "output tree not initialised");\r
+    return 11;\r
+  }\r
+\r
+  fParam->SetMaxTBin(fRecoParam->GetLastBin());//set number of timebins from reco -> param\r
+  //   used later for memory allocation\r
+\r
+  Bool_t isAltro=kFALSE;\r
+\r
+  AliTPCROC * roc = AliTPCROC::Instance();\r
+  AliTPCCalPad * noiseTPC = AliTPCcalibDB::Instance()->GetPadNoise();\r
+  AliTPCAltroMapping** mapping =AliTPCcalibDB::Instance()->GetMapping();\r
+  //\r
+  AliTPCRawStream input(rawReader,(AliAltroMapping**)mapping);\r
+\r
+  const Int_t kNIS = fParam->GetNInnerSector();//number of inner sectors\r
+  const Int_t kNOS = fParam->GetNOuterSector();//number of outer sectors\r
+  const Int_t kNS = kNIS + kNOS;//all sectors\r
+\r
+\r
+  //crate TPC view\r
+  AliTPCDigitsArray *digarr=new AliTPCDigitsArray(kFALSE);//data not sim\r
+  digarr->Setup(fParam);//as usually parameters\r
+\r
+  //\r
+  // Loop over sectors\r
+  //\r
+  for(Int_t iSec = 0; iSec < kNS; iSec++) {\r
+    AliTPCCalROC * noiseROC;\r
+    if(noiseTPC==0x0){\r
+      noiseROC =new AliTPCCalROC(iSec);//noise=0\r
+    }\r
+    else{\r
+      noiseROC = noiseTPC->GetCalROC(iSec);  // noise per given sector\r
+    }\r
+    Int_t nRows = 0; //number of rows in sector\r
+    Int_t nDDLs = 0; //number of DDLs\r
+    Int_t indexDDL = 0; //DDL index\r
+    if (iSec < kNIS) {\r
+      nRows = fParam->GetNRowLow();\r
+      nDDLs = 2;\r
+      indexDDL = iSec * 2;\r
+    }else {\r
+      nRows = fParam->GetNRowUp();\r
+      nDDLs = 4;\r
+      indexDDL = (iSec-kNIS) * 4 + kNIS * 2;\r
+    }\r
+\r
+    //\r
+    // Load the raw data for corresponding DDLs\r
+    //\r
+    rawReader->Reset();\r
+    rawReader->Select("TPC",indexDDL,indexDDL+nDDLs-1);\r
+\r
+    if(input.Next()) {\r
+      isAltro=kTRUE;\r
+      // Allocate memory for rows in sector (pads(depends on row) x timebins)\r
+      for(Int_t iRow = 0; iRow < nRows; iRow++) {\r
+       digarr->CreateRow(iSec,iRow); \r
+      }//end loop over rows\r
+    }\r
+    rawReader->Reset();\r
+    rawReader->Select("TPC",indexDDL,indexDDL+nDDLs-1);\r
+\r
+    //\r
+    // Begin loop over altro data\r
+    //\r
+    while (input.Next()) {\r
+      \r
+      //check sector consistency\r
+      if (input.GetSector() != iSec)\r
+       AliFatal(Form("Sector index mismatch ! Expected (%d), but got (%d) !",iSec,input.GetSector()));\r
+      \r
+      Int_t iRow = input.GetRow();\r
+      Int_t iPad = input.GetPad();\r
+      Int_t iTimeBin = input.GetTime();\r
+\r
+      //\r
+      if(fDebugLevel==72){\r
+       fHistoRow->Fill(iRow);\r
+       fHistoPad->Fill(iPad);\r
+       fHistoTime->Fill(iTimeBin);\r
+       fHistoRowPad->Fill(iPad,iRow);\r
+      }else if(fDebugLevel>=0&&fDebugLevel<72){\r
+       if(iSec==fDebugLevel){\r
+         fHistoRow->Fill(iRow);\r
+         fHistoPad->Fill(iPad);\r
+         fHistoTime->Fill(iTimeBin);\r
+         fHistoRowPad->Fill(iPad,iRow);\r
+       }\r
+      }else if(fDebugLevel==73){\r
+       if(iSec<36){\r
+         fHistoRow->Fill(iRow);\r
+         fHistoPad->Fill(iPad);\r
+         fHistoTime->Fill(iTimeBin);\r
+         fHistoRowPad->Fill(iPad,iRow);\r
+       }\r
+      }else if(fDebugLevel==74){\r
+       if(iSec>=36){\r
+         fHistoRow->Fill(iRow);\r
+         fHistoPad->Fill(iPad);\r
+         fHistoTime->Fill(iTimeBin);\r
+         fHistoRowPad->Fill(iPad,iRow);\r
+       }\r
+      }\r
+\r
+      //check row consistency\r
+      if (iRow < 0 ) continue;\r
+      if (iRow < 0 || iRow >= nRows){\r
+       AliError(Form("Pad-row index (%d) outside the range (%d -> %d) !",\r
+                     iRow, 0, nRows -1));\r
+       continue;\r
+      }\r
+\r
+      //check pad consistency\r
+      if (iPad < 0 || iPad >= (Int_t)(roc->GetNPads(iSec,iRow))) {\r
+       AliError(Form("Pad index (%d) outside the range (%d -> %d) !",\r
+                     iPad, 0, roc->GetNPads(iSec,iRow) ));\r
+       continue;\r
+      }\r
+\r
+      //check time consistency\r
+      if ( iTimeBin < fRecoParam->GetFirstBin() || iTimeBin >= fRecoParam->GetLastBin()){\r
+       //cout<<iTimeBin<<endl;\r
+       continue;\r
+       AliFatal(Form("Timebin index (%d) outside the range (%d -> %d) !",\r
+                     iTimeBin, 0, fRecoParam->GetLastBin() -1));\r
+      }\r
+\r
+      //signal\r
+      Int_t signal = input.GetSignal();\r
+      if (signal <= fZeroSup || \r
+         iTimeBin < fFirstBin ||\r
+         iTimeBin > fLastBin\r
+         ) {\r
+       digarr->GetRow(iSec,iRow)->SetDigitFast(0,iTimeBin,iPad);\r
+       continue;\r
+      }\r
+      if (!noiseROC) continue;\r
+      Double_t noiseOnPad = noiseROC->GetValue(iRow,iPad);//noise on given pad and row in sector\r
+      if (noiseOnPad > fMaxNoiseAbs){\r
+       digarr->GetRow(iSec,iRow)->SetDigitFast(0,iTimeBin,iPad);\r
+       continue; // consider noisy pad as dead\r
+      }\r
+      if(signal <= fMaxNoiseSigma * noiseOnPad){\r
+       digarr->GetRow(iSec,iRow)->SetDigitFast(0,iTimeBin,iPad);\r
+       continue;\r
+      }\r
+      digarr->GetRow(iSec,iRow)->SetDigitFast(signal,iTimeBin,iPad);      \r
+    }//end of loop over altro data\r
+  }//end of loop over sectors\r
+  \r
+  SetDigArr(digarr);\r
+  if(isAltro) FindClusterKrIO();\r
+  delete digarr;\r
+\r
+  return 0;\r
+}\r
+\r
+void AliTPCclustererKr::CleanSector(Int_t sector){\r
+  //\r
+  // clean isolated digits\r
+  //  \r
+  const Int_t kNRows=fParam->GetNRow(sector);//number of rows in sector\r
+  for(Int_t iRow=0; iRow<kNRows; ++iRow){\r
+    AliSimDigits *digrow;\r
+    if(fRawData){\r
+      digrow = (AliSimDigits*)fDigarr->GetRow(sector,iRow);//real data\r
+    }else{\r
+      digrow = (AliSimDigits*)fDigarr->LoadRow(sector,iRow);//MC\r
+    }\r
+    if(!digrow) continue;\r
+    digrow->ExpandBuffer(); //decrunch\r
+    const Int_t kNPads = digrow->GetNCols();  // number of pads\r
+    const Int_t kNTime = digrow->GetNRows(); // number of timebins\r
+    for(Int_t iPad=1;iPad<kNPads-1;iPad++){\r
+      Short_t*  val = digrow->GetDigitsColumn(iPad);\r
+\r
+      for(Int_t iTimeBin=1;iTimeBin<kNTime-1;iTimeBin++){\r
+       if (val[iTimeBin]<=0) continue;\r
+       if (val[iTimeBin-1]+val[iTimeBin+1]<fIsolCut) {val[iTimeBin]=0; continue;}\r
+       if (val[iTimeBin-kNTime]+val[iTimeBin+kNTime]<fIsolCut) {val[iTimeBin]=0; continue;}\r
+       //\r
+       if (val[iTimeBin-1-kNTime]+val[iTimeBin+1+kNTime]<fIsolCut) {val[iTimeBin]=0; continue;}\r
+       if (val[iTimeBin+1-kNTime]+val[iTimeBin-1+kNTime]<fIsolCut) {val[iTimeBin]=0; continue;}\r
+\r
+      }\r
+    }\r
+  }\r
+}\r
+\r
+\r
+////____________________________________________________________________________\r
+Int_t AliTPCclustererKr::FindClusterKrIO()\r
+{\r
+\r
+  //\r
+  //fParam and  fDigarr must be set to run this method\r
+  //\r
+\r
+  Int_t clusterCounter=0;\r
+  const Int_t nTotalSector=fParam->GetNSector();//number of sectors\r
+  for(Int_t iSec=0; iSec<nTotalSector; ++iSec){\r
+    CleanSector(iSec);\r
+\r
+    //vector of maxima for each sector\r
+    //std::vector<AliPadMax*> maximaInSector;\r
+    TObjArray *maximaInSector=new TObjArray();//to store AliPadMax*\r
+\r
+    //\r
+    //  looking for the maxima on the pad\r
+    //\r
+\r
+    const Int_t kNRows=fParam->GetNRow(iSec);//number of rows in sector\r
+    for(Int_t iRow=0; iRow<kNRows; ++iRow){\r
+      AliSimDigits *digrow;\r
+      if(fRawData){\r
+       digrow = (AliSimDigits*)fDigarr->GetRow(iSec,iRow);//real data\r
+      }else{\r
+       digrow = (AliSimDigits*)fDigarr->LoadRow(iSec,iRow);//MC\r
+      }\r
+      if(digrow){//if pointer exist\r
+       digrow->ExpandBuffer(); //decrunch\r
+       const Int_t kNPads = digrow->GetNCols();  // number of pads\r
+       const Int_t kNTime = digrow->GetNRows(); // number of timebins\r
+       for(Int_t iPad=0;iPad<kNPads;iPad++){\r
+         \r
+         Int_t timeBinMax=-1;//timebin of maximum \r
+         Int_t valueMaximum=-1;//value of maximum in adc\r
+         Int_t increaseBegin=-1;//timebin when increase starts\r
+         Int_t sumAdc=0;//sum of adc on the pad in maximum surrounding\r
+         bool ifIncreaseBegin=true;//flag - check if increasing started\r
+         bool ifMaximum=false;//flag - check if it could be maximum\r
+         Short_t* val = digrow->GetDigitsColumn(iPad);\r
+         for(Int_t iTimeBin=1;iTimeBin<kNTime-1;iTimeBin++){\r
+           if (!ifMaximum)  {\r
+             if (val[iTimeBin]==-1) break;   // 0 until the end\r
+             for( ; iTimeBin<kNTime-2&&val[iTimeBin]<fMinAdc ;iTimeBin++);\r
+           }\r
+           //\r
+           Short_t adc = val[iTimeBin];\r
+\r
+           if(adc<fMinAdc){//standard was 3 for fMinAdc\r
+             if(ifMaximum){\r
+               if(iTimeBin-increaseBegin<fMinTimeBins){//at least 2 time bins\r
+                 timeBinMax=-1;\r
+                 valueMaximum=-1;\r
+                 increaseBegin=-1;\r
+                 sumAdc=0;\r
+                 ifIncreaseBegin=true;\r
+                 ifMaximum=false;\r
+                 continue;\r
+               }\r
+               //insert maximum, default values and set flags\r
+               //Double_t xCord,yCord;\r
+               //GetXY(iSec,iRow,iPad,xCord,yCord);\r
+               Double_t x[]={iRow,iPad,iTimeBin};\r
+               Int_t i[]={iSec};\r
+               AliTPCTransform trafo;\r
+               trafo.Transform(x,i,0,1);\r
+               \r
+               AliPadMax *oneMaximum = new AliPadMax(AliTPCvtpr(valueMaximum,\r
+                                                                timeBinMax,\r
+                                                                iPad,\r
+                                                                iRow,\r
+                                                                x[0],//xCord,\r
+                                                                x[1],//yCord,\r
+                                                                x[2]/*timeBinMax*/),\r
+                                                     increaseBegin,\r
+                                                     iTimeBin-1,\r
+                                                     sumAdc);\r
+               maximaInSector->AddLast(oneMaximum);\r
+               \r
+               timeBinMax=-1;\r
+               valueMaximum=-1;\r
+               increaseBegin=-1;\r
+               sumAdc=0;\r
+               ifIncreaseBegin=true;\r
+               ifMaximum=false;\r
+             }\r
+             continue;\r
+           }\r
+\r
+\r
+\r
+\r
+\r
+\r
+           if(ifIncreaseBegin){\r
+             ifIncreaseBegin=false;\r
+             increaseBegin=iTimeBin;\r
+           }\r
+           \r
+           if(adc>valueMaximum){\r
+             timeBinMax=iTimeBin;\r
+             valueMaximum=adc;\r
+             ifMaximum=true;\r
+           }\r
+           sumAdc+=adc;\r
+           if(iTimeBin==kNTime-1 && ifMaximum && kNTime-increaseBegin>fMinTimeBins){//on the edge\r
+             //at least 3 timebins\r
+             //insert maximum, default values and set flags\r
+             //Double_t xCord,yCord;\r
+             //GetXY(iSec,iRow,iPad,xCord,yCord);\r
+             Double_t x[]={iRow,iPad,iTimeBin};\r
+             Int_t i[]={iSec};\r
+             AliTPCTransform trafo;\r
+             trafo.Transform(x,i,0,1);\r
+             AliPadMax *oneMaximum = new AliPadMax(AliTPCvtpr(valueMaximum,\r
+                                                              timeBinMax,\r
+                                                              iPad,\r
+                                                              iRow,\r
+                                                              x[0],//xCord,\r
+                                                              x[1],//yCord,\r
+                                                              x[2]/*timeBinMax*/),\r
+                                                   increaseBegin,\r
+                                                   iTimeBin-1,\r
+                                                   sumAdc);\r
+             maximaInSector->AddLast(oneMaximum);\r
+               \r
+             timeBinMax=-1;\r
+             valueMaximum=-1;\r
+             increaseBegin=-1;\r
+             sumAdc=0;\r
+             ifIncreaseBegin=true;\r
+             ifMaximum=false;\r
+             continue;\r
+           }\r
+           \r
+         }//end loop over timebins\r
+       }//end loop over pads\r
+//      }else{\r
+//     cout<<"Pointer does not exist!!"<<endl;\r
+      }//end if poiner exists\r
+    }//end loop over rows\r
+\r
+    MakeClusters(maximaInSector,iSec,clusterCounter);\r
+    //\r
+    maximaInSector->SetOwner(kTRUE);\r
+    maximaInSector->Delete();\r
+    delete maximaInSector;\r
+  }//end sector for\r
+  cout<<"Number of clusters in event: "<<clusterCounter<<endl;\r
+  return 0;\r
+}\r
+\r
+void AliTPCclustererKr::MakeClusters(TObjArray * maximaInSector, Int_t iSec, Int_t &clusterCounter){\r
+  //\r
+  // Make clusters\r
+  //\r
+\r
+  Int_t maxDig=0;\r
+  Int_t maxSumAdc=0;\r
+  Int_t maxTimeBin=0;\r
+  Int_t maxPad=0;\r
+  Int_t maxRow=0;\r
+  Double_t maxX=0;\r
+  Double_t maxY=0;\r
+  Double_t maxT=0;\r
+  Int_t entriesArr = maximaInSector->GetEntriesFast();\r
+  for(Int_t it1 = 0; it1 < entriesArr; ++it1 ) {\r
+    \r
+    AliPadMax *mp1=(AliPadMax *)maximaInSector->UncheckedAt(it1);\r
+    if (!mp1) continue;\r
+    AliTPCclusterKr clusterKr;\r
+    \r
+    Int_t nUsedPads=1;\r
+    Int_t clusterValue=0;\r
+    clusterValue+=(mp1)->GetSum();\r
+    list<Int_t> nUsedRows;\r
+    nUsedRows.push_back((mp1)->GetRow());\r
+    \r
+    maxDig      =(mp1)->GetAdc() ;\r
+    maxSumAdc   =(mp1)->GetSum() ;\r
+    maxTimeBin  =(mp1)->GetTime();\r
+    maxPad      =(mp1)->GetPad() ;\r
+    maxRow      =(mp1)->GetRow() ;\r
+    maxX        =(mp1)->GetX();\r
+    maxY        =(mp1)->GetY();\r
+    maxT        =(mp1)->GetT();\r
+    \r
+    AliSimDigits *digrowTmp;\r
+    if(fRawData){\r
+      digrowTmp = (AliSimDigits*)fDigarr->GetRow(iSec,(mp1)->GetRow());\r
+    }else{\r
+      digrowTmp = (AliSimDigits*)fDigarr->LoadRow(iSec,(mp1)->GetRow());\r
+    }\r
+    \r
+    digrowTmp->ExpandBuffer(); //decrunch\r
+    \r
+    for(Int_t itb=(mp1)->GetBegin(); itb<((mp1)->GetEnd())+1; itb++){\r
+      Int_t adcTmp = digrowTmp->GetDigitUnchecked(itb,(mp1)->GetPad());\r
+      AliTPCvtpr *vtpr=new AliTPCvtpr(adcTmp,itb,(mp1)->GetPad(),(mp1)->GetRow(),(mp1)->GetX(),(mp1)->GetY(),(mp1)->GetT());\r
+      clusterKr.AddDigitToCluster(vtpr);\r
+    }\r
+    clusterKr.SetCenter();//set centr of the cluster\r
+    \r
+    for(Int_t it2 = it1+1; it2 < entriesArr; ++it2 ) {\r
+      AliPadMax *mp2=(AliPadMax *)maximaInSector->UncheckedAt(it2);\r
+      if (!mp2) continue;\r
+      if (TMath::Abs(clusterKr.GetCenterX() - (mp2)->GetX()) > fMaxPadRangeCm) continue;\r
+      if (TMath::Abs(clusterKr.GetCenterY() - (mp2)->GetY()) > fMaxRowRangeCm) continue;      \r
+      if (TMath::Abs(clusterKr.GetCenterT() - (mp2)->GetT()) > fMaxTimeRange) continue;\r
+\r
+      {\r
+       clusterValue+=(mp2)->GetSum();\r
+       \r
+       nUsedPads++;\r
+       nUsedRows.push_back((mp2)->GetRow());\r
+       \r
+       AliSimDigits *digrowTmp1;\r
+       if(fRawData){\r
+         digrowTmp1 = (AliSimDigits*)fDigarr->GetRow(iSec,(mp2)->GetRow());\r
+       }else{\r
+         digrowTmp1 = (AliSimDigits*)fDigarr->LoadRow(iSec,(mp2)->GetRow());\r
+       }\r
+       digrowTmp1->ExpandBuffer(); //decrunch\r
+       \r
+       for(Int_t itb=(mp2)->GetBegin(); itb<(mp2)->GetEnd()+1; itb++){\r
+         Int_t adcTmp = digrowTmp1->GetDigitUnchecked(itb,(mp2)->GetPad());\r
+         AliTPCvtpr *vtpr=new AliTPCvtpr(adcTmp,itb,(mp2)->GetPad(),(mp2)->GetRow(),(mp2)->GetX(),(mp2)->GetY(),(mp2)->GetT());\r
+         clusterKr.AddDigitToCluster(vtpr);\r
+       }\r
+       \r
+       clusterKr.SetCenter();//set center of the cluster\r
+       \r
+       //which one is bigger\r
+       if( (mp2)->GetAdc() > maxDig ){\r
+         maxDig      =(mp2)->GetAdc() ;\r
+         maxSumAdc   =(mp2)->GetSum() ;\r
+         maxTimeBin  =(mp2)->GetTime();\r
+         maxPad      =(mp2)->GetPad() ;\r
+         maxRow      =(mp2)->GetRow() ;\r
+         maxX        =(mp2)->GetX() ;\r
+         maxY        =(mp2)->GetY() ;\r
+         maxT        =(mp2)->GetT() ;\r
+       } else if ( (mp2)->GetAdc() == maxDig ){\r
+         if( (mp2)->GetSum() > maxSumAdc){\r
+           maxDig      =(mp2)->GetAdc() ;\r
+           maxSumAdc   =(mp2)->GetSum() ;\r
+           maxTimeBin  =(mp2)->GetTime();\r
+           maxPad      =(mp2)->GetPad() ;\r
+           maxRow      =(mp2)->GetRow() ;\r
+           maxX        =(mp2)->GetX() ;\r
+           maxY        =(mp2)->GetY() ;\r
+           maxT        =(mp2)->GetT() ;\r
+         }\r
+       }\r
+       delete maximaInSector->RemoveAt(it2);\r
+      }\r
+    }//inside loop\r
+    delete maximaInSector->RemoveAt(it1);          \r
+    clusterKr.SetSize();\r
+    //through out clusters on the edge and noise\r
+    //if(clusterValue/clusterKr.fCluster.size()<fValueToSize)continue;\r
+    if(clusterValue/(clusterKr.GetSize())<fValueToSize)continue;\r
+    \r
+    clusterKr.SetADCcluster(clusterValue);\r
+    clusterKr.SetNPads(nUsedPads);\r
+    clusterKr.SetMax(AliTPCvtpr(maxDig,maxTimeBin,maxPad,maxRow,maxX,maxY,maxT));\r
+    clusterKr.SetSec(iSec);\r
+    clusterKr.SetSize();\r
+    \r
+    nUsedRows.sort();\r
+    nUsedRows.unique();\r
+    clusterKr.SetNRows(nUsedRows.size());\r
+    clusterKr.SetCenter();\r
+    \r
+    clusterKr.SetRMS();//Set pad,row,timebin RMS\r
+    clusterKr.Set1D();//Set size in pads and timebins\r
+\r
+    clusterCounter++;\r
+    \r
+    \r
+    //save each cluster into file\r
+    if (fOutput){\r
+      (*fOutput)<<"Kr"<<\r
+       "Cl.="<<&clusterKr<<\r
+       "\n";\r
+    }\r
+    //end of save each cluster into file adc.root\r
+  }//outer loop\r
+}\r
+\r
+\r
+\r
+////____________________________________________________________________________\r
+\r
+\r
+void AliTPCclustererKr::GetXY(Int_t sec,Int_t row,Int_t pad,Double_t& xGlob,Double_t& yGlob){\r
+  //\r
+  //gives global XY coordinate of the pad\r
+  //\r
+\r
+  Double_t yLocal = fParam->GetPadRowRadii(sec,row);//radius of row in sector in cm\r
+\r
+  Int_t padmax = fParam->GetNPads(sec,row);//number of pads in a given row\r
+  Float_t padXSize;\r
+  if(sec<fParam->GetNInnerSector())padXSize=0.4;\r
+  else padXSize=0.6;\r
+  Double_t xLocal=(pad+0.5-padmax/2.)*padXSize;//x-value of the center of pad\r
+\r
+  Float_t sin,cos;\r
+  fParam->AdjustCosSin((Int_t)sec,cos,sin);//return sinus and cosinus of the sector\r
+\r
+  Double_t xGlob1 =  xLocal * cos + yLocal * sin;\r
+  Double_t yGlob1 = -xLocal * sin + yLocal * cos;\r
+\r
+\r
+  Double_t rot=0;\r
+  rot=TMath::Pi()/2.;\r
+\r
+  xGlob =  xGlob1 * TMath::Cos(rot) + yGlob1 * TMath::Sin(rot);\r
+  yGlob = -xGlob1 * TMath::Sin(rot) + yGlob1 * TMath::Cos(rot);\r
+\r
+   yGlob=-1*yGlob;\r
+   if(sec<18||(sec>=36&&sec<54)) xGlob =-1*xGlob;\r
+\r
+\r
+  return;\r
+}\r
index 86fce7809b03a715094691718bad15568a9e7145..46a01d10e0429ac210434555145d7fd6805de91b 100644 (file)
-#ifndef ALITPCCLUSTERERKR_H
-#define ALITPCCLUSTERERKR_H
-/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
- * See cxx source for full Copyright notice                               */
-
-/* $Id: AliTPCclusterKr.h,v 1.8 2008/02/07 16:07:15 matyja Exp $ */
-
-//-------------------------------------------------------
-//                    TPC Kr Clusterer Class
-//
-//   Origin: Adam Matyja, INP PAN, adam.matyja@ifj.edu.pl
-//-------------------------------------------------------
-
-#include "TObject.h"
-
-class AliTPCclusterKr;
-class AliPadMax;
-class AliSimDigits;
-class AliTPCv4;
-class AliTPCParam;
-class AliTPCDigitsArray;
-class AliTPCvtpr;
-class AliTPCClustersRow;
-class TTree;
-class TH1F;
-class TH2F;
-
-class AliTPCTransform;
-
-//used in raw data finder
-class AliTPCROC;
-class AliTPCCalPad;
-class AliTPCAltroMapping;
-class AliTPCcalibDB;
-class AliTPCRawStream;
-class AliTPCRecoParam;
-class AliTPCReconstructor;
-class AliRawReader;
-class AliTPCCalROC;
-class TTreeSRedirector;
-//_____________________________________________________________________________
-class AliTPCclustererKr: public TObject{
-public:
-  AliTPCclustererKr();
-  AliTPCclustererKr(const AliTPCclustererKr &param);//copy constructor
-  AliTPCclustererKr &operator = (const AliTPCclustererKr & param); 
-  virtual ~AliTPCclustererKr();
-
-  //finders
-  virtual Int_t FinderIO();//for MC
-  virtual Int_t FinderIO(AliRawReader* rawReader);//for data
-  virtual Int_t FindClusterKrIO();//main routine for finding clusters
-  virtual Int_t CleanSector(Int_t sector); // clean isolated digits
-
-
-  //other
-  void GetXY(Int_t sec,Int_t row,Int_t pad,Double_t& xGlob,Double_t& yGlob);//give XY coordinate of the pad
-
-  virtual void SetInput(TTree * tree);  //set input tree with digits    
-  virtual void SetOutput(TTree * tree); //set output tree with clusters
-  virtual void SetParam(AliTPCParam *param){fParam=param;}//set TPC parameters
-  virtual void SetDigArr(AliTPCDigitsArray *digarr){fDigarr=digarr;}//set current array of digits
-  virtual void SetRecoParam(AliTPCRecoParam *recoParam=0);//set reconstruction parameters
-
-
-
-  //setters for cluster finder parameters
-  virtual void SetZeroSup(Int_t v){fZeroSup=v;}//set zero suppresion parameter
-  virtual void SetFirstBin(Int_t v){fFirstBin=v;}//set first considered timebin
-  virtual void SetLastBin(Int_t v){fLastBin=v;}//set last considered timebin
-  virtual void SetMaxNoiseAbs(Float_t v){fMaxNoiseAbs=v;}//set maximal noise value
-  virtual void SetMaxNoiseSigma(Float_t v){fMaxNoiseSigma=v;}//set maximal noise sigma
-
-  virtual void SetMinAdc(Int_t v){v<=0?fMinAdc=1:fMinAdc=v;}//set fMinAdc
-  virtual void SetMinTimeBins(Int_t v){fMinTimeBins=v;}//set fMinTimeBins
-//  virtual void SetMaxPadRange(Int_t v){fMaxPadRange=v;}//set fMaxPadRange
-//  virtual void SetMaxRowRange(Int_t v){fMaxRowRange=v;}//set fMaxRowRange
-  virtual void SetMaxTimeRange(Int_t v){fMaxTimeRange=v;}//set fMaxTimeRange
-  virtual void SetValueToSize(Float_t v){fValueToSize=v;}//set fValueToSize
-
-  virtual void SetMaxPadRangeCm(Double_t v){fMaxPadRangeCm=v;}//set fMaxPadRangeCm
-  virtual void SetMaxRowRangeCm(Double_t v){fMaxRowRangeCm=v;}//set fMaxRowRangeCm
-
-  virtual void SetDebugLevel(Int_t debug){fDebugLevel=debug;}
-  //debug = 0 to 71 -sector number to  print
-  // = 72 - all sectors
-  // = 73 - inners
-  // = 74 - outers
-
-  virtual void SetHistoRow(TH1F *histo)   {fHistoRow   =histo;}
-  virtual void SetHistoPad(TH1F *histo)   {fHistoPad   =histo;}
-  virtual void SetHistoTime(TH1F *histo)  {fHistoTime  =histo;}
-  virtual void SetHistoRowPad(TH2F *histo){fHistoRowPad=histo;}
-
-  //getters for cluster finder parameters
-  Int_t GetZeroSup() const {return fZeroSup;}//get zero suppresion parameter
-  Int_t GetFirstBin() const {return fFirstBin;}//get first considered timebin
-  Int_t GetLastBin() const {return fLastBin;}//get last considered timebin
-  Float_t GetMaxNoiseAbs() const {return fMaxNoiseAbs;}//get maximal noise value
-  Float_t GetMaxNoiseSigma() const {return fMaxNoiseSigma;}//get maximal noise sigma
-
-  Int_t GetMinAdc() const {return fMinAdc;}//get fMinAdc
-  Int_t GetMinTimeBins() const {return fMinTimeBins;}//get fMinTimeBins
-//  Int_t GetMaxPadRange() const {return fMaxPadRange;}//get fMaxPadRange
-//  Int_t GetMaxRowRange() const {return fMaxRowRange;}//get fMaxRowRange
-  Int_t GetMaxTimeRange() const {return fMaxTimeRange;}//get fMaxTimeRange
-  Float_t GetValueToSize() const {return fValueToSize;}//get fValueToSize
-
-  Double_t GetMaxPadRangeCm() const {return fMaxPadRangeCm;}//get fMaxPadRangeCm
-  Double_t GetMaxRowRangeCm() const {return fMaxRowRangeCm;}//get fMaxRowRangeCm
-
-  Int_t GetDebugLevel() const {return fDebugLevel;}
-  TH1F * GetHistoRow(){return fHistoRow;}
-  TH1F * GetHistoPad(){return fHistoPad;}
-  TH1F * GetHistoTime(){return fHistoTime;}
-  TH2F * GetHistoRowPad(){return fHistoRowPad;}
-
-private:
-  void MakeClusters(TObjArray * maximaInSector, Int_t iSec, Int_t &clusterCounter);
-  Bool_t fRawData; //flag =0 for MC =1 for real data
-
-  TTree * fInput;   //!input  tree with digits - object not owner
-  TTreeSRedirector * fOutput;   //!output tree with clusters - object not owner
-
-  AliTPCParam * fParam;//!TPC parameters
-  AliTPCDigitsArray *fDigarr;//! pointer to current array if digits
-
-  //only for raw data :)
-  const AliTPCRecoParam  * fRecoParam;        //! reconstruction parameters
-
-  //cluster finder parameters
-  Int_t fZeroSup;//zero suppresion parameter = 2 def.
-  Int_t fFirstBin;//first considered time bin used by cluster finder = 60 def.
-  Int_t fLastBin;//last considered time bin used by cluster finder = 950 def.
-  Float_t fMaxNoiseAbs;// maximal noise value on pad used in cluster finder = 2 def.
-  Float_t fMaxNoiseSigma;// maximal noise sigma on pad used in cluster finder = 3 def.
-
-  Int_t fMinAdc;//minimal value of acd count in each timebin = 3 def.
-  Int_t fMinTimeBins;//minimal value of time bins one after each other = 2 def.
-//  Int_t fMaxPadRange;//maximal pad range from maximum = 4 def.
-//  Int_t fMaxRowRange;//maximal row range from maximum = 3 def.
-  Int_t fMaxTimeRange;//maximal time bin range from maximum = 7 def.
-  Float_t fValueToSize;//ratio cluster value to cluster size = 3.5 def.
-
-  Double_t fMaxPadRangeCm;//maximal pad range in cm from maximum = 2.5cm def.
-  Double_t fMaxRowRangeCm;//maximal row range in cm from maximum = 3.5cm def.
-
-  Int_t fDebugLevel;//! debug level variable
-  TH1F *fHistoRow;//!debug histo for rows
-  TH1F *fHistoPad;//!debug histo for pads
-  TH1F *fHistoTime;//!debug histo for timebins
-  TH2F *fHistoRowPad;//!debug histo for rows and pads
-
-  ClassDef(AliTPCclustererKr,6)  // Time Projection Chamber Kr clusters
-};
-
-
-#endif
-
-
+#ifndef ALITPCCLUSTERERKR_H\r
+#define ALITPCCLUSTERERKR_H\r
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *\r
+ * See cxx source for full Copyright notice                               */\r
+\r
+/* $Id: AliTPCclusterKr.h,v 1.8 2008/02/07 16:07:15 matyja Exp $ */\r
+\r
+//-------------------------------------------------------\r
+//                    TPC Kr Clusterer Class\r
+//\r
+//   Origin: Adam Matyja, INP PAN, adam.matyja@ifj.edu.pl\r
+//-------------------------------------------------------\r
+\r
+#include "TObject.h"\r
+\r
+class AliTPCclusterKr;\r
+class AliPadMax;\r
+class AliSimDigits;\r
+class AliTPCv4;\r
+class AliTPCParam;\r
+class AliTPCDigitsArray;\r
+class AliTPCvtpr;\r
+class AliTPCClustersRow;\r
+class TTree;\r
+class TH1F;\r
+class TH2F;\r
+\r
+class AliTPCTransform;\r
+\r
+//used in raw data finder\r
+class AliTPCROC;\r
+class AliTPCCalPad;\r
+class AliTPCAltroMapping;\r
+class AliTPCcalibDB;\r
+class AliTPCRawStream;\r
+class AliTPCRecoParam;\r
+class AliTPCReconstructor;\r
+class AliRawReader;\r
+class AliTPCCalROC;\r
+class TTreeSRedirector;\r
+//_____________________________________________________________________________\r
+class AliTPCclustererKr: public TObject{\r
+public:\r
+  AliTPCclustererKr();\r
+  AliTPCclustererKr(const AliTPCclustererKr &param);//copy constructor\r
+  AliTPCclustererKr &operator = (const AliTPCclustererKr & param); \r
+  virtual ~AliTPCclustererKr();\r
+\r
+  //finders\r
+  virtual Int_t FinderIO();//for MC\r
+  virtual Int_t FinderIO(AliRawReader* rawReader);//for data\r
+  virtual Int_t FindClusterKrIO();//main routine for finding clusters\r
+  virtual void CleanSector(Int_t sector); // clean isolated digits\r
+\r
+  //other\r
+  void GetXY(Int_t sec,Int_t row,Int_t pad,Double_t& xGlob,Double_t& yGlob);//give XY coordinate of the pad\r
+\r
+  virtual void SetInput(TTree * tree);  //set input tree with digits    \r
+  virtual void SetOutput(TTree * tree); //set output tree with clusters\r
+  virtual void SetParam(AliTPCParam *param){fParam=param;}//set TPC parameters\r
+  virtual void SetDigArr(AliTPCDigitsArray *digarr){fDigarr=digarr;}//set current array of digits\r
+  virtual void SetRecoParam(AliTPCRecoParam *recoParam=0);//set reconstruction parameters\r
+\r
+\r
+\r
+  //setters for cluster finder parameters\r
+  virtual void SetZeroSup(Int_t v){fZeroSup=v;}//set zero suppresion parameter\r
+  virtual void SetFirstBin(Int_t v){fFirstBin=v;}//set first considered timebin\r
+  virtual void SetLastBin(Int_t v){fLastBin=v;}//set last considered timebin\r
+  virtual void SetMaxNoiseAbs(Float_t v){fMaxNoiseAbs=v;}//set maximal noise value\r
+  virtual void SetMaxNoiseSigma(Float_t v){fMaxNoiseSigma=v;}//set maximal noise sigma\r
+\r
+  virtual void SetMinAdc(Int_t v){v<=0?fMinAdc=1:fMinAdc=v;}//set fMinAdc\r
+  virtual void SetMinTimeBins(Int_t v){fMinTimeBins=v;}//set fMinTimeBins\r
+//  virtual void SetMaxPadRange(Int_t v){fMaxPadRange=v;}//set fMaxPadRange\r
+//  virtual void SetMaxRowRange(Int_t v){fMaxRowRange=v;}//set fMaxRowRange\r
+  virtual void SetMaxTimeRange(Int_t v){fMaxTimeRange=v;}//set fMaxTimeRange\r
+  virtual void SetValueToSize(Float_t v){fValueToSize=v;}//set fValueToSize\r
+\r
+  virtual void SetMaxPadRangeCm(Double_t v){fMaxPadRangeCm=v;}//set fMaxPadRangeCm\r
+  virtual void SetMaxRowRangeCm(Double_t v){fMaxRowRangeCm=v;}//set fMaxRowRangeCm\r
+\r
+  virtual void SetIsolCut(Short_t v){fIsolCut=v;}\r
+\r
+  virtual void SetDebugLevel(Int_t debug){fDebugLevel=debug;}\r
+  //debug = 0 to 71 -sector number to  print\r
+  // = 72 - all sectors\r
+  // = 73 - inners\r
+  // = 74 - outers\r
+\r
+  virtual void SetHistoRow(TH1F *histo)   {fHistoRow   =histo;}\r
+  virtual void SetHistoPad(TH1F *histo)   {fHistoPad   =histo;}\r
+  virtual void SetHistoTime(TH1F *histo)  {fHistoTime  =histo;}\r
+  virtual void SetHistoRowPad(TH2F *histo){fHistoRowPad=histo;}\r
+\r
+  //getters for cluster finder parameters\r
+  Int_t GetZeroSup() const {return fZeroSup;}//get zero suppresion parameter\r
+  Int_t GetFirstBin() const {return fFirstBin;}//get first considered timebin\r
+  Int_t GetLastBin() const {return fLastBin;}//get last considered timebin\r
+  Float_t GetMaxNoiseAbs() const {return fMaxNoiseAbs;}//get maximal noise value\r
+  Float_t GetMaxNoiseSigma() const {return fMaxNoiseSigma;}//get maximal noise sigma\r
+\r
+  Int_t GetMinAdc() const {return fMinAdc;}//get fMinAdc\r
+  Int_t GetMinTimeBins() const {return fMinTimeBins;}//get fMinTimeBins\r
+//  Int_t GetMaxPadRange() const {return fMaxPadRange;}//get fMaxPadRange\r
+//  Int_t GetMaxRowRange() const {return fMaxRowRange;}//get fMaxRowRange\r
+  Int_t GetMaxTimeRange() const {return fMaxTimeRange;}//get fMaxTimeRange\r
+  Float_t GetValueToSize() const {return fValueToSize;}//get fValueToSize\r
+\r
+  Double_t GetMaxPadRangeCm() const {return fMaxPadRangeCm;}//get fMaxPadRangeCm\r
+  Double_t GetMaxRowRangeCm() const {return fMaxRowRangeCm;}//get fMaxRowRangeCm\r
+\r
+  Short_t GetIsolCut() const {return fIsolCut;}\r
+\r
+  Int_t GetDebugLevel() const {return fDebugLevel;}\r
+  TH1F * GetHistoRow(){return fHistoRow;} \r
+  TH1F * GetHistoPad(){return fHistoPad;}\r
+  TH1F * GetHistoTime(){return fHistoTime;}\r
+  TH2F * GetHistoRowPad(){return fHistoRowPad;}\r
+\r
+private:\r
+  void MakeClusters(TObjArray * maximaInSector, Int_t iSec, Int_t &clusterCounter);\r
+  Bool_t fRawData; //flag =0 for MC =1 for real data\r
+\r
+  TTree * fInput;   //!input  tree with digits - object not owner\r
+  TTreeSRedirector * fOutput;   //!output tree with clusters - object not owner\r
+\r
+  AliTPCParam * fParam;//!TPC parameters\r
+  AliTPCDigitsArray *fDigarr;//! pointer to current array if digits\r
+\r
+  //only for raw data :)\r
+  const AliTPCRecoParam  * fRecoParam;        //! reconstruction parameters\r
+\r
+  //cluster finder parameters\r
+  Int_t fZeroSup;//zero suppresion parameter = 2 def.\r
+  Int_t fFirstBin;//first considered time bin used by cluster finder = 60 def.\r
+  Int_t fLastBin;//last considered time bin used by cluster finder = 950 def.\r
+  Float_t fMaxNoiseAbs;// maximal noise value on pad used in cluster finder = 2 def.\r
+  Float_t fMaxNoiseSigma;// maximal noise sigma on pad used in cluster finder = 3 def.\r
+\r
+  Int_t fMinAdc;//minimal value of acd count in each timebin = 3 def.\r
+  Int_t fMinTimeBins;//minimal value of time bins one after each other = 2 def.\r
+//  Int_t fMaxPadRange;//maximal pad range from maximum = 4 def.\r
+//  Int_t fMaxRowRange;//maximal row range from maximum = 3 def.\r
+  Int_t fMaxTimeRange;//maximal time bin range from maximum = 7 def.\r
+  Float_t fValueToSize;//ratio cluster value to cluster size = 3.5 def.\r
+\r
+  Double_t fMaxPadRangeCm;//maximal pad range in cm from maximum = 2.5cm def.\r
+  Double_t fMaxRowRangeCm;//maximal row range in cm from maximum = 3.5cm def.\r
+\r
+  Short_t fIsolCut;//isolation cut in 3D = 5 def.\r
+\r
+  Int_t fDebugLevel;//! debug level variable\r
+  TH1F *fHistoRow;//!debug histo for rows\r
+  TH1F *fHistoPad;//!debug histo for pads\r
+  TH1F *fHistoTime;//!debug histo for timebins\r
+  TH2F *fHistoRowPad;//!debug histo for rows and pads\r
+\r
+  ClassDef(AliTPCclustererKr,6)  // Time Projection Chamber Kr clusters\r
+};\r
+\r
+\r
+#endif\r
+\r
+\r
index d68596967e689a9afa689e8c2261bb0310caf6a5..d00b9ccc07f7b79713194aca6348566bbe5ea036 100644 (file)
-/****************************************************************************
- *           Origin: A.Matyja amatyja@cern.ch                               *
- ****************************************************************************/
-
-/*
-
-  macro to create array of clusters from TPC digits
-  input files - galice.root 
-                digits.root - file with digits - usualy use link to galice.root
-                           - in splitted mode - neccesary to create link to proper file
-                           
-   output file - TPC.RecPoints.root
-
-
-//  Warning - if cluster file AliTPCclusters.root already exist - macro exit and don't produce anything
-              
-*/
-
-
-#ifndef __CINT__
-#include <iostream.h>
-#include "AliRun.h"
-#include "AliTPCv4.h"
-#include "AliTPCParam.h"
-#include "AliTPCclusterKr.h"
-#include "AliTPCclustererKr.h"
-#include "TFile.h"
-#include "TStopwatch.h"
-#include "TTree.h"
-#endif
-
-Int_t FindKrClusters(){
-
-  //
-  //Load DataBase
-  //
-  char *ocdbpath ="local:///afs/cern.ch/alice/tpctest/OCDB";
-  //char *ocdbpath ="local:///home/matyja/baza/OCDB";
-  if (ocdbpath==0){
-    ocdbpath="alien://folder=/alice/data/2007/LHC07w/OCDB/";
-  }
-  printf("OCDB PATH = %s\n",ocdbpath); 
-  AliCDBManager * man = AliCDBManager::Instance();
-  man->SetDefaultStorage(ocdbpath);
-  man->SetRun(0);
-
-  AliRunLoader* rl = AliRunLoader::Open("galice.root");
-  if (rl == 0x0) {
-    cerr<<"Can not open session"<<endl;
-    return 1;
-  }
-  
-  AliTPCLoader *tpcl = (AliTPCLoader*)rl->GetLoader("TPCLoader");
-  if (tpcl == 0x0) {
-    cerr<<"Can not get TPC Loader"<<endl;
-    return 1;
-  }
-
-  if (tpcl->LoadDigits()) {
-    cerr<<"Error occured while loading digits"<<endl;
-    return 1;
-  }
-
-  if (rl->LoadgAlice()) {
-    cerr<<"Error occured while LoadgAlice"<<endl;
-    return 1;
-  }
-  
-  gAlice=rl->GetAliRun();
-  if (!gAlice) {
-    cerr<<"Can't get gAlice !\n";
-    return 1;
-  }
-
-  TDirectory *cwd = gDirectory;
-
-  AliTPCv4 *tpc = (AliTPCv4*)gAlice->GetDetector("TPC");
-  Int_t ver = tpc->IsVersion(); 
-  cerr<<"TPC version "<<ver<<" has been found !\n";
-
-  rl->CdGAFile();
-  
-  AliTPCParam *param=(AliTPCParamSR *)gDirectory->Get("75x40_100x60_150x60");
-  if (!param) {cerr<<"TPC parameters have not been found !\n"; return 4;}
-  
-  AliTPCDigitsArray *digarr=new AliTPCDigitsArray;
-  digarr->Setup(param);
-
-  cerr<<"It has begun"<<endl;  
-  TStopwatch timer;
-  timer.Start();
-  cwd->cd();
-  
-  Int_t nevmax=rl->GetNumberOfEvents();//number of events in run
-  for(Int_t nev=0;nev<nevmax /*&& nev<1*/ ;nev++){
-    rl->GetEvent(nev);
-    
-    TTree* input_tree= tpcl->TreeD();//tree with digits
-    if (input_tree == 0x0){
-      cerr << "Can not get TreeD for event " <<nev<<endl;
-      continue;
-    }
-    
-    digarr->ConnectTree(input_tree);
-    
-    TTree *output_tree =tpcl->TreeR();
-    if(output_tree==0x0){
-      tpcl->MakeTree("R");
-      output_tree = tpcl->TreeR();
-      if (output_tree == 0x0){
-       cerr << "Problems with output tree (TreeR) for event "<<nev<<endl;
-       continue;
-      }
-    }
-    
-    cout<<"Processing event "<<nev<<endl;
-
-    //test
-    //cout<<"nentries"<<input_tree->GetEntries();
-
-    AliTPCclustererKr *clusters = new AliTPCclustererKr();
-    clusters->SetParam(param);
-    clusters->SetInput(input_tree);
-    clusters->SetOutput(output_tree);
-    clusters->SetDigArr(digarr);
-    clusters->FinderIO();
-
-    tpcl->WriteRecPoints("OVERWRITE");
-  }
-  
-
-  timer.Stop(); timer.Print();
-  
-  delete rl;//cleans everything
-  
-  return 0;
-}
+\r
+/****************************************************************************\r
+ *           Origin: A.Matyja amatyja@cern.ch                               *\r
+ ****************************************************************************/\r
+\r
+/*\r
+\r
+  macro to create array of clusters from TPC digits\r
+  input files - galice.root \r
+                digits.root - file with digits - usualy use link to galice.root\r
+                           - in splitted mode - neccesary to create link to proper file\r
+                           \r
+   output file - TPC.RecPoints.root\r
+\r
+\r
+//  Warning - if cluster file AliTPCclusters.root already exist - macro exit and don't produce anything\r
+              \r
\r
+*/\r
+\r
+\r
+#ifndef __CINT__\r
+#include <iostream.h>\r
+#include "AliRun.h"\r
+#include "AliTPCv4.h"\r
+#include "AliTPCParam.h"\r
+#include "AliTPCclusterKr.h"\r
+#include "AliTPCclustererKr.h"\r
+#include "TFile.h"\r
+#include "TStopwatch.h"\r
+#include "TTree.h"\r
+#endif\r
+\r
+Int_t FindKrClusters(){\r
+\r
+  //\r
+  //Load DataBase\r
+  //\r
+  //char *ocdbpath ="local:///afs/cern.ch/alice/tpctest/OCDB";\r
+  //char *ocdbpath ="local:///home/matyja/baza/OCDB";\r
+  char *ocdbpath ="local:///data/baza/OCDB";\r
+  if (ocdbpath==0){\r
+    ocdbpath="alien://folder=/alice/data/2007/LHC07w/OCDB/";\r
+  }\r
+  printf("OCDB PATH = %s\n",ocdbpath); \r
+  AliCDBManager * man = AliCDBManager::Instance();\r
+  man->SetDefaultStorage(ocdbpath);\r
+  man->SetRun(0);\r
+\r
+  AliRunLoader* rl = AliRunLoader::Open("galice.root");\r
+  if (rl == 0x0) {\r
+    cerr<<"Can not open session"<<endl;\r
+    return 1;\r
+  }\r
+  \r
+  AliTPCLoader *tpcl = (AliTPCLoader*)rl->GetLoader("TPCLoader");\r
+  if (tpcl == 0x0) {\r
+    cerr<<"Can not get TPC Loader"<<endl;\r
+    return 1;\r
+  }\r
+\r
+  if (tpcl->LoadDigits()) {\r
+    cerr<<"Error occured while loading digits"<<endl;\r
+    return 1;\r
+  }\r
+\r
+  if (rl->LoadgAlice()) {\r
+    cerr<<"Error occured while LoadgAlice"<<endl;\r
+    return 1;\r
+  }\r
+  \r
+  gAlice=rl->GetAliRun();\r
+  if (!gAlice) {\r
+    cerr<<"Can't get gAlice !\n";\r
+    return 1;\r
+  }\r
+\r
+  TDirectory *cwd = gDirectory;\r
+\r
+  AliTPCv4 *tpc = (AliTPCv4*)gAlice->GetDetector("TPC");\r
+  Int_t ver = tpc->IsVersion(); \r
+  cerr<<"TPC version "<<ver<<" has been found !\n";\r
+\r
+  rl->CdGAFile();\r
+  \r
+  AliTPCParam *param=(AliTPCParamSR *)gDirectory->Get("75x40_100x60_150x60");\r
+  if (!param) {cerr<<"TPC parameters have not been found !\n"; return 4;}\r
+  \r
+  AliTPCDigitsArray *digarr=new AliTPCDigitsArray;\r
+  digarr->Setup(param);\r
+\r
+  cerr<<"It has begun"<<endl;  \r
+  TStopwatch timer;\r
+  timer.Start();\r
+  cwd->cd();\r
+\r
+  TTree *output_tree;\r
+  \r
+  AliTPCclustererKr *clusters = new AliTPCclustererKr();\r
+  clusters->SetParam(param);\r
+  clusters->SetOutput(output_tree);\r
+\r
+  clusters->SetMinAdc(3);//signal threshold (everything below is treated as 0)\r
+  clusters->SetMinTimeBins(2);//number of neighbouring timebins\r
+  clusters->SetMaxPadRangeCm(5.);//distance of the cluster center to the center of a pad (in cm)\r
+  clusters->SetMaxRowRangeCm(5.);//distance of the cluster center to the center of a padrow (in cm)\r
+  clusters->SetMaxTimeRange(7.);//distance of the cluster center to the max time bin on a pad (in tackts)\r
+  //ie. fabs(centerT - time)<7\r
+\r
+  clusters->SetIsolCut(3);//set isolation cut threshold\r
+  clusters->SetValueToSize(3.1);//cut reduce peak at 0\r
+\r
+  Int_t nevmax=rl->GetNumberOfEvents();//number of events in run\r
+  for(Int_t nev=0;nev<nevmax /*&& nev<1*/ ;nev++){\r
+    rl->GetEvent(nev);\r
+    \r
+    TTree* input_tree= tpcl->TreeD();//tree with digits\r
+    if (input_tree == 0x0){\r
+      cerr << "Can not get TreeD for event " <<nev<<endl;\r
+      continue;\r
+    }\r
+    digarr->ConnectTree(input_tree);\r
+    clusters->SetInput(input_tree);\r
+    clusters->SetDigArr(digarr);\r
+    cout<<"Processing event "<<nev<<endl;\r
+    clusters->FinderIO();\r
+\r
+  }\r
+  delete clusters;\r
+  timer.Stop(); timer.Print();\r
+  \r
+  delete rl;//cleans everything\r
+  \r
+  return 0;\r
+}\r
index 0569bd9b7033d1f92ebca7bb37afa3683558c766..0a00aeadc249fdd91527b513b50f8cf05e389c34 100644 (file)
-//
-
-Int_t FindKrClusterCheck(const char *fileName="data.root");
-
-
-Int_t FindKrClustersRaw(const char *fileName="data.root"){
-
-  
-
-  //
-  // remove Altro warnings
-  //
-  AliLog::SetClassDebugLevel("AliTPCRawStream",-5);
-  AliLog::SetClassDebugLevel("AliRawReaderDate",-5);
-  AliLog::SetClassDebugLevel("AliTPCAltroMapping",-5);
-  AliLog::SetModuleDebugLevel("RAW",-5);
-  //
-  // Get calibration
-  //
-  char *ocdbpath = gSystem->Getenv("OCDB_PATH");
-  //char *ocdbpath ="local:///afs/cern.ch/alice/tpctest/OCDB";
-  if (ocdbpath==0){
-    ocdbpath="alien://folder=/alice/data/2007/LHC07w/OCDB/";
-  }
-  printf("OCDB PATH = %s\n",ocdbpath); 
-  AliCDBManager * man = AliCDBManager::Instance();
-  man->SetDefaultStorage(ocdbpath);
-  man->SetRun(100000000);
-
-  AliTPCCalPad * noiseTPC = AliTPCcalibDB::Instance()->GetPadNoise();
-  AliTPCAltroMapping** mapping =AliTPCcalibDB::Instance()->GetMapping();
-  //
-  //
-
-
-  //define tree
-  TFile *hfile=new TFile("adc.root","RECREATE","ADC file");
-  // Create a ROOT Tree
-  TTree *mytree = new TTree("Kr","Krypton cluster tree");
-
-
-  Int_t debugLevel=1;
-  if(debugLevel>0){
-    TH1F *histoRow   =new TH1F("histoRow","rows",100,0.,100.);
-    TH1F *histoPad   =new TH1F("histoPad","pads",150,0.,150.);
-    TH1F *histoTime  =new TH1F("histoTime","timebins",100,0.,1000.);
-    TH2F *histoRowPad=new TH2F("histoRowPad","pads-vs-rows" ,150,0.,150.,100,0.,100.);
-  }
-
-
-  //one general output
-  AliTPCclustererKr *clusters = new AliTPCclustererKr();
-  clusters->SetOutput(mytree);
-  clusters->SetRecoParam(0);
-
-  if(debugLevel>0){
-    clusters->SetDebugLevel(debugLevel);
-    clusters->SetHistoRow(histoRow   );
-    clusters->SetHistoPad(histoPad   );
-    clusters->SetHistoTime(histoTime  );
-    clusters->SetHistoRowPad(histoRowPad);
-  }
-
-
-  AliTPCParamSR *param=new AliTPCParamSR();
-  //only for geometry parameters loading - temporarly
-//  AliRunLoader* rl = AliRunLoader::Open("galice.root");
-//  AliTPCParam *param=(AliTPCParamSR *)gDirectory->Get("75x40_100x60_150x60");
-  //if (!param) {cerr<<"TPC parameters have not been found !\n"; return 4;}
-
-  clusters->SetParam(param);
-
-  //set cluster finder parameters (from data)
-  clusters->SetZeroSup(param->GetZeroSup());//zero suppression parameter
-  clusters->SetFirstBin(60);//first bin
-  clusters->SetLastBin(950);//last bin
-  clusters->SetMaxNoiseAbs(2);//maximal noise
-  clusters->SetMaxNoiseSigma(3);//maximal amount of sigma of noise
-
-  //set cluster finder parameters (from MC)
-  clusters->SetMinAdc(3);//signal threshold (everything below is treated as 0)
-  clusters->SetMinTimeBins(2);//number of neighbouring timebins
-  clusters->SetMaxPadRangeCm(2.5);//distance of the cluster center to the center of a pad (in cm)
-  clusters->SetMaxRowRangeCm(3.5);//distance of the cluster center to the center of a padrow (in cm)
-  clusters->SetMaxTimeRange(7);//distance of the cluster center to the max time bin on a pad (in tackts)
-  //ie. fabs(centerT - time)<7
-  clusters->SetValueToSize(6);//cut reduce peak at 0
-
-
-  AliRawReader *reader = new AliRawReaderRoot(fileName);
-  reader->Reset();
-
-  TStopwatch timer;
-  timer.Start();
-
-  AliAltroRawStreamFast* stream = new AliAltroRawStreamFast(reader);
-  stream->SelectRawData("TPC");
-
-  Int_t evtnr=0;
-  while (reader->NextEvent()) {
-    //output for each event
-
-    //if(evtnr>4) break;
-    cout<<"Evt = "<<evtnr<<endl;
-    clusters->FinderIO(reader);
-    evtnr++;
-    AliSysInfo::AddStamp(Form("Event%d",evtnr),evtnr);
-  }
-
-
-  mytree->Print();//print rootuple summary 
-  // Save all objects in this file
-  hfile->Write();
-  // Close the file
-  hfile->Close();
-
-  timer.Stop();
-  timer.Print();
-  printf("Deleting clusterer\n");
-  delete clusters;
-  printf("Deleting stream\n");
-  delete stream;
-  printf("Deleting raw reader\n");
-  delete reader;
-
-//   TCanvas *c2=new TCanvas("c2","title",800,800);
-//   c2->SetHighLightColor(2);
-//   c2->Range(-458.9552,-2948.238,3296.642,26856.6);
-//   c2->SetBorderSize(2);
-//   c2->SetLeftMargin(0.15);
-//   c2->SetRightMargin(0.06);
-//   c2->SetFrameFillColor(0);
-
-//   gStyle->SetOptStat(111111);
-//   histoRow->Draw();
-//   c2->Print("rows.ps");
-//   histoPad->Draw();
-//   c2->Print("pads.ps");
-//   histoTime->Draw();
-//   c2->Print("timebins.ps");
-//   histoRowPad->Draw();
-//   c2->Print("row-pad.ps");
-
-  return 0;
-}
-
-
-Int_t FindKrClusterCheck(const char *fileName){
-  //
-  //
-  gSystem->Load("$ROOTSYS/lib/libGui.so");
-  gSystem->Load("$ROOTSYS/lib/libTree.so");
-  gSystem->Load("$MEMSTAT/libMemStat.so");
-  {
-    TMemStat memstat(1000000,100000,kTRUE);
-    AliSysInfo::AddCallBack(TMemStatManager::GetInstance()->fStampCallBack);
-    FindKrClustersRaw(fileName); 
-  }
-  // the output memstat.root file
-  TMemStat draw("memstat.root");
-  // Print some information
-  // code info
-  draw.MakeReport(0,0);
-
-}
+//\r
+\r
+Int_t FindKrClusterCheck(const char *fileName="data.root");\r
+\r
+\r
+Int_t FindKrClustersRaw(const char *fileName="data.root"){\r
+\r
+  \r
+\r
+  //\r
+  // remove Altro warnings\r
+  //\r
+  AliLog::SetClassDebugLevel("AliTPCRawStream",-5);\r
+  AliLog::SetClassDebugLevel("AliRawReaderDate",-5);\r
+  AliLog::SetClassDebugLevel("AliTPCAltroMapping",-5);\r
+  AliLog::SetModuleDebugLevel("RAW",-5);\r
+  //\r
+  // Get calibration\r
+  //\r
+  char *ocdbpath = gSystem->Getenv("OCDB_PATH");\r
+  //char *ocdbpath ="local:///afs/cern.ch/alice/tpctest/OCDB";\r
+  if (ocdbpath==0){\r
+    ocdbpath="alien://folder=/alice/data/2007/LHC07w/OCDB/";\r
+  }\r
+  printf("OCDB PATH = %s\n",ocdbpath); \r
+  AliCDBManager * man = AliCDBManager::Instance();\r
+  man->SetDefaultStorage(ocdbpath);\r
+  man->SetRun(100000000);\r
+\r
+  AliTPCCalPad * noiseTPC = AliTPCcalibDB::Instance()->GetPadNoise();\r
+  AliTPCAltroMapping** mapping =AliTPCcalibDB::Instance()->GetMapping();\r
+  //\r
+  //\r
+\r
+\r
+  //define tree\r
+  TFile *hfile=new TFile("adc.root","RECREATE","ADC file");\r
+  // Create a ROOT Tree\r
+  TTree *mytree = new TTree("Kr","Krypton cluster tree");\r
+\r
+\r
+  Int_t debugLevel=1;\r
+  if(debugLevel>0){\r
+    TH1F *histoRow   =new TH1F("histoRow","rows",100,0.,100.);\r
+    TH1F *histoPad   =new TH1F("histoPad","pads",150,0.,150.);\r
+    TH1F *histoTime  =new TH1F("histoTime","timebins",100,0.,1000.);\r
+    TH2F *histoRowPad=new TH2F("histoRowPad","pads-vs-rows" ,150,0.,150.,100,0.,100.);\r
+  }\r
+\r
+\r
+  //one general output\r
+  AliTPCclustererKr *clusters = new AliTPCclustererKr();\r
+  clusters->SetOutput(mytree);\r
+  clusters->SetRecoParam(0);\r
+\r
+  if(debugLevel>0){\r
+    clusters->SetDebugLevel(debugLevel);\r
+    clusters->SetHistoRow(histoRow   );\r
+    clusters->SetHistoPad(histoPad   );\r
+    clusters->SetHistoTime(histoTime  );\r
+    clusters->SetHistoRowPad(histoRowPad);\r
+  }\r
+\r
+\r
+  AliTPCParamSR *param=new AliTPCParamSR();\r
+  //only for geometry parameters loading - temporarly\r
+//  AliRunLoader* rl = AliRunLoader::Open("galice.root");\r
+//  AliTPCParam *param=(AliTPCParamSR *)gDirectory->Get("75x40_100x60_150x60");\r
+  //if (!param) {cerr<<"TPC parameters have not been found !\n"; return 4;}\r
+\r
+  clusters->SetParam(param);\r
+\r
+  //set cluster finder parameters (from data)\r
+  clusters->SetZeroSup(param->GetZeroSup());//zero suppression parameter\r
+  clusters->SetFirstBin(60);//first bin\r
+  clusters->SetLastBin(950);//last bin\r
+  clusters->SetMaxNoiseAbs(2);//maximal noise\r
+  clusters->SetMaxNoiseSigma(3);//maximal amount of sigma of noise\r
+\r
+  //set cluster finder parameters (from MC)\r
+  clusters->SetMinAdc(3);//signal threshold (everything below is treated as 0)\r
+  clusters->SetMinTimeBins(2);//number of neighbouring timebins\r
+  clusters->SetMaxPadRangeCm(5.);//distance of the cluster center to the center of a pad (in cm)\r
+  clusters->SetMaxRowRangeCm(5.);//distance of the cluster center to the center of a padrow (in cm)\r
+  clusters->SetMaxTimeRange(5.);//distance of the cluster center to the max time bin on a pad (in tackts)\r
+  //ie. fabs(centerT - time)<7\r
+  clusters->SetValueToSize(7.);//cut reduce peak at 0\r
+  clusters->SetIsolCut(3);//set isolation cut threshold\r
+\r
+  AliRawReader *reader = new AliRawReaderRoot(fileName);\r
+  reader->Reset();\r
+\r
+  TStopwatch timer;\r
+  timer.Start();\r
+\r
+  AliAltroRawStreamFast* stream = new AliAltroRawStreamFast(reader);\r
+  stream->SelectRawData("TPC");\r
+\r
+  Int_t evtnr=0;\r
+  while (reader->NextEvent()) {\r
+    //output for each event\r
+\r
+    //if(evtnr>4) break;\r
+    cout<<"Evt = "<<evtnr<<endl;\r
+    clusters->FinderIO(reader);\r
+    evtnr++;\r
+    AliSysInfo::AddStamp(Form("Event%d",evtnr),evtnr);\r
+  }\r
+\r
+\r
+  mytree->Print();//print rootuple summary \r
+  // Save all objects in this file\r
+  hfile->Write();\r
+  // Close the file\r
+  hfile->Close();\r
+\r
+  timer.Stop();\r
+  timer.Print();\r
+  printf("Deleting clusterer\n");\r
+  delete clusters;\r
+  printf("Deleting stream\n");\r
+  delete stream;\r
+  printf("Deleting raw reader\n");\r
+  delete reader;\r
+\r
+//   TCanvas *c2=new TCanvas("c2","title",800,800);\r
+//   c2->SetHighLightColor(2);\r
+//   c2->Range(-458.9552,-2948.238,3296.642,26856.6);\r
+//   c2->SetBorderSize(2);\r
+//   c2->SetLeftMargin(0.15);\r
+//   c2->SetRightMargin(0.06);\r
+//   c2->SetFrameFillColor(0);\r
+\r
+//   gStyle->SetOptStat(111111);\r
+//   histoRow->Draw();\r
+//   c2->Print("rows.ps");\r
+//   histoPad->Draw();\r
+//   c2->Print("pads.ps");\r
+//   histoTime->Draw();\r
+//   c2->Print("timebins.ps");\r
+//   histoRowPad->Draw();\r
+//   c2->Print("row-pad.ps");\r
+\r
+  return 0;\r
+}\r
+\r
+\r
+Int_t FindKrClusterCheck(const char *fileName){\r
+  //\r
+  //\r
+  gSystem->Load("$ROOTSYS/lib/libGui.so");\r
+  gSystem->Load("$ROOTSYS/lib/libTree.so");\r
+  gSystem->Load("$MEMSTAT/libMemStat.so");\r
+  {\r
+    TMemStat memstat(1000000,100000,kTRUE);\r
+    AliSysInfo::AddCallBack(TMemStatManager::GetInstance()->fStampCallBack);\r
+    FindKrClustersRaw(fileName); \r
+  }\r
+  // the output memstat.root file\r
+  TMemStat draw("memstat.root");\r
+  // Print some information\r
+  // code info\r
+  draw.MakeReport(0,0);\r
+\r
+}\r