-
/**************************************************************************
* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
* *
#include "AliTRDdataArrayF.h"
#include "AliTRDdataArrayI.h"
#include "AliTRDdataArrayS.h"
+#include "AliTRDdataArrayDigits.h"
#include "AliTRDdigitsManager.h"
#include "AliTRDrawData.h"
#include "AliTRDcalibDB.h"
#include "Cal/AliTRDCalROC.h"
#include "Cal/AliTRDCalDet.h"
+#include "Cal/AliTRDCalSingleChamberStatus.h"
ClassImp(AliTRDclusterizer)
//_____________________________________________________________________________
-AliTRDclusterizer::AliTRDclusterizer()
+AliTRDclusterizer::AliTRDclusterizer(AliTRDReconstructor *rec)
:TNamed()
+ ,fReconstructor(rec)
,fRunLoader(NULL)
,fClusterTree(NULL)
,fRecPoints(NULL)
+ ,fTrackletTree(NULL)
,fDigitsManager(NULL)
+ ,fTrackletContainer(NULL)
,fAddLabels(kTRUE)
,fRawVersion(2)
,fIndexesOut(NULL)
// AliTRDclusterizer default constructor
//
+ AliTRDcalibDB *trd = 0x0;
+ if (!(trd = AliTRDcalibDB::Instance())) {
+ AliFatal("Could not get calibration object");
+ }
+
fRawVersion = AliTRDfeeParam::Instance()->GetRAWversion();
+ // Initialize debug stream
+ if(fReconstructor){
+ if(fReconstructor->GetStreamLevel(AliTRDReconstructor::kClusterizer) > 1){
+ TDirectory *savedir = gDirectory;
+ //fgDebugStreamer = new TTreeSRedirector("TRD.ClusterizerDebug.root");
+ savedir->cd();
+ }
+ }
}
//_____________________________________________________________________________
-AliTRDclusterizer::AliTRDclusterizer(const Text_t *name, const Text_t *title)
+AliTRDclusterizer::AliTRDclusterizer(const Text_t *name, const Text_t *title, AliTRDReconstructor *rec)
:TNamed(name,title)
+ ,fReconstructor(rec)
,fRunLoader(NULL)
,fClusterTree(NULL)
,fRecPoints(NULL)
+ ,fTrackletTree(NULL)
,fDigitsManager(new AliTRDdigitsManager())
+ ,fTrackletContainer(NULL)
,fAddLabels(kTRUE)
,fRawVersion(2)
,fIndexesOut(NULL)
// AliTRDclusterizer constructor
//
+ AliTRDcalibDB *trd = 0x0;
+ if (!(trd = AliTRDcalibDB::Instance())) {
+ AliFatal("Could not get calibration object");
+ }
+
+ // Initialize debug stream
+ if(fReconstructor){
+ if(fReconstructor->GetStreamLevel(AliTRDReconstructor::kClusterizer) > 1){
+ TDirectory *savedir = gDirectory;
+ //fgDebugStreamer = new TTreeSRedirector("TRD.ClusterizerDebug.root");
+ savedir->cd();
+ }
+ }
+
fDigitsManager->CreateArrays();
fRawVersion = AliTRDfeeParam::Instance()->GetRAWversion();
//_____________________________________________________________________________
AliTRDclusterizer::AliTRDclusterizer(const AliTRDclusterizer &c)
:TNamed(c)
+ ,fReconstructor(c.fReconstructor)
,fRunLoader(NULL)
,fClusterTree(NULL)
,fRecPoints(NULL)
+ ,fTrackletTree(NULL)
+ ,fTrackletContainer(NULL)
,fDigitsManager(NULL)
,fAddLabels(kTRUE)
,fRawVersion(2)
fDigitsManager = NULL;
}
+ if (fTrackletContainer)
+ {
+ delete fTrackletContainer;
+ fTrackletContainer = NULL;
+ }
+
if (fIndexesOut)
{
delete fIndexesOut;
{
((AliTRDclusterizer &) c).Copy(*this);
}
+
return *this;
}
((AliTRDclusterizer &) c).fClusterTree = NULL;
((AliTRDclusterizer &) c).fRecPoints = NULL;
+ ((AliTRDclusterizer &) c).fTrackletTree = NULL;
((AliTRDclusterizer &) c).fDigitsManager = NULL;
+ ((AliTRDclusterizer &) c).fTrackletContainer = NULL;
((AliTRDclusterizer &) c).fAddLabels = fAddLabels;
((AliTRDclusterizer &) c).fRawVersion = fRawVersion;
((AliTRDclusterizer &) c).fIndexesOut = NULL;
fClusterTree = clusterTree;
fClusterTree->Branch("TRDcluster","TObjArray",&ioArray,32000,0);
+
+ // tracklet writing
+ if (fReconstructor->IsWritingTracklets()){
+ TString evfoldname = AliConfig::GetDefaultEventFolderName();
+ fRunLoader = AliRunLoader::GetRunLoader(evfoldname);
+
+ if (!fRunLoader) {
+ fRunLoader = AliRunLoader::Open("galice.root");
+ }
+ if (!fRunLoader) {
+ AliError(Form("Can not open session for file galice.root."));
+ return kFALSE;
+ }
+
+ UInt_t **leaves = new UInt_t *[2];
+ AliDataLoader *dl = fRunLoader->GetLoader("TRDLoader")->GetDataLoader("tracklets");
+ if (!dl) {
+ AliError("Could not get the tracklets data loader!");
+ AliDataLoader *dl = new AliDataLoader("TRD.Tracklets.root","tracklets", "tracklets");
+ fRunLoader->GetLoader("TRDLoader")->AddDataLoader(dl);
+ }
+ else {
+ fTrackletTree = dl->Tree();
+ if (!fTrackletTree)
+ {
+ dl->MakeTree();
+ fTrackletTree = dl->Tree();
+ }
+ TBranch *trkbranch = fTrackletTree->GetBranch("trkbranch");
+ if (!trkbranch)
+ fTrackletTree->Branch("trkbranch",leaves[0],"det/i:side/i:tracklets[256]/i");
+ }
+ }
+
return kTRUE;
}
}
+//_____________________________________________________________________________
+Bool_t AliTRDclusterizer::WriteTracklets(Int_t det)
+{
+
+ UInt_t **leaves = new UInt_t *[2];
+ for (Int_t i=0; i<2 ;i++){
+ leaves[i] = new UInt_t[258];
+ leaves[i][0] = det; // det
+ leaves[i][1] = i; // side
+ memcpy(leaves[i]+2, fTrackletContainer[i], sizeof(UInt_t) * 256);
+ }
+
+
+ if (!fTrackletTree){
+ AliDataLoader *dl = fRunLoader->GetLoader("TRDLoader")->GetDataLoader("tracklets");
+ dl->MakeTree();
+ fTrackletTree = dl->Tree();
+ }
+
+ TBranch *trkbranch = fTrackletTree->GetBranch("trkbranch");
+ if (!trkbranch) {
+ trkbranch = fTrackletTree->Branch("trkbranch",leaves[0],"det/i:side/i:tracklets[256]/i");
+ }
+
+ for (Int_t i=0; i<2; i++){
+ if (leaves[i][2]>0) {
+ trkbranch->SetAddress(leaves[i]);
+ fTrackletTree->Fill();
+ }
+ }
+
+ AliDataLoader *dl = fRunLoader->GetLoader("TRDLoader")->GetDataLoader("tracklets");
+ dl->WriteData("OVERWRITE");
+ //dl->Unload();
+ delete [] leaves;
+
+ return kTRUE;
+}
+
//_____________________________________________________________________________
void AliTRDclusterizer::ResetHelperIndexes(AliTRDSignalIndex *indexesIn)
{
for (Int_t i = 0; i < AliTRDgeometry::kNdet; i++)
{
- AliTRDdataArrayS *digitsIn = (AliTRDdataArrayS *) fDigitsManager->GetDigits(i);
+ AliTRDdataArrayDigits *digitsIn = (AliTRDdataArrayDigits*) fDigitsManager->GetDigits(i);
// This is to take care of switched off super modules
if (!digitsIn->HasData())
{
fDigitsManager->SetUseDictionaries(fAddLabels);
+ // tracklet container for raw tracklet writing
+ if (!fTrackletContainer && fReconstructor->IsWritingTracklets())
+ {
+ fTrackletContainer = new UInt_t *[2];
+ for (Int_t i=0; i<2 ;i++){
+ fTrackletContainer[i] = new UInt_t[256]; // maximum tracklets for one HC
+ }
+ }
+
AliTRDrawStreamBase *pinput = AliTRDrawStreamBase::GetRawStream(rawReader);
AliTRDrawStreamBase &input = *pinput;
AliInfo(Form("Stream version: %s", input.IsA()->GetName()));
Int_t det = 0;
- while ((det = input.NextChamber(fDigitsManager)) >= 0)
+ while ((det = input.NextChamber(fDigitsManager,fTrackletContainer)) >= 0)
{
Bool_t iclusterBranch = kFALSE;
if (fDigitsManager->GetIndexes(det)->HasEntry())
- {
- iclusterBranch = MakeClusters(det);
- }
+ {
+ iclusterBranch = MakeClusters(det);
+ }
if (iclusterBranch == kFALSE)
- {
- WriteClusters(det);
- ResetRecPoints();
- }
+ {
+ WriteClusters(det);
+ ResetRecPoints();
+ }
fDigitsManager->RemoveDigits(det);
- fDigitsManager->RemoveDictionaries(det);
+ fDigitsManager->RemoveDictionaries(det);
fDigitsManager->ClearIndexes(det);
+
+ if (!fReconstructor->IsWritingTracklets()) continue;
+ if (*(fTrackletContainer[0]) > 0 || *(fTrackletContainer[1]) > 0) WriteTracklets(det); // if there is tracklet words in this det
}
+
+ if (fReconstructor->IsWritingTracklets()){
+ delete [] fTrackletContainer;
+ fTrackletContainer = NULL;
+ }
delete fDigitsManager;
fDigitsManager = NULL;
}
+//_____________________________________________________________________________
+UChar_t AliTRDclusterizer::GetStatus(Short_t &signal)
+{
+ //
+ // Check if a pad is masked
+ //
+
+ UChar_t status = 0;
+
+ if(signal>0 && TESTBIT(signal, 10)){
+ CLRBIT(signal, 10);
+ for(int ibit=0; ibit<4; ibit++){
+ if(TESTBIT(signal, 11+ibit)){
+ SETBIT(status, ibit);
+ CLRBIT(signal, 11+ibit);
+ }
+ }
+ }
+ return status;
+}
+
//_____________________________________________________________________________
Bool_t AliTRDclusterizer::MakeClusters(Int_t det)
{
// Get the digits
// digits should be expanded beforehand!
// digitsIn->Expand();
- AliTRDdataArrayS *digitsIn = (AliTRDdataArrayS *) fDigitsManager->GetDigits(det);
-
+ AliTRDdataArrayDigits *digitsIn = (AliTRDdataArrayDigits *) fDigitsManager->GetDigits(det);
+
// This is to take care of switched off super modules
if (!digitsIn->HasData())
{
// There is no ADC threshold anymore, and simParam should not be used in clusterizer. KO
Float_t adcThreshold = 0;
+ if (!fReconstructor){
+ AliError("Reconstructor not set\n");
+ return kFALSE;
+ }
+
// Threshold value for the maximum
- Float_t maxThresh = AliTRDReconstructor::RecoParam()->GetClusMaxThresh();
+ Float_t maxThresh = fReconstructor->GetRecoParam() ->GetClusMaxThresh();
// Threshold value for the digit signal
- Float_t sigThresh = AliTRDReconstructor::RecoParam()->GetClusSigThresh();
+ Float_t sigThresh = fReconstructor->GetRecoParam() ->GetClusSigThresh();
+
+ // Threshold value for the maximum ( cut noise)
+ Float_t minMaxCutSigma = fReconstructor->GetRecoParam() ->GetMinMaxCutSigma();
+ // Threshold value for the sum pad ( cut noise)
+ Float_t minLeftRightCutSigma = fReconstructor->GetRecoParam() ->GetMinLeftRightCutSigma();
// Iteration limit for unfolding procedure
const Float_t kEpsilon = 0.01;
Double_t padSignal[kNsig];
Double_t clusterSignal[kNclus];
- Int_t icham = indexesIn->GetChamber();
- Int_t iplan = indexesIn->GetPlane();
- Int_t isect = indexesIn->GetSM();
+ Int_t istack = indexesIn->GetStack();
+ Int_t ilayer = indexesIn->GetLayer();
+ Int_t isector = indexesIn->GetSM();
// Start clustering in the chamber
- Int_t idet = AliTRDgeometry::GetDetector(iplan,icham,isect);
+ Int_t idet = AliTRDgeometry::GetDetector(ilayer,istack,isector);
if (idet != det)
{
AliError("Strange Detector number Missmatch!");
// TRD space point transformation
fTransform->SetDetector(det);
- Int_t ilayer = AliGeomManager::kTRD1 + iplan;
- Int_t imodule = icham + AliTRDgeometry::Ncham() * isect;
- UShort_t volid = AliGeomManager::LayerToVolUID(ilayer,imodule);
+ Int_t iGeoLayer = AliGeomManager::kTRD1 + ilayer;
+ Int_t iGeoModule = istack + AliTRDgeometry::Nstack() * isector;
+ UShort_t volid = AliGeomManager::LayerToVolUID(iGeoLayer,iGeoModule);
Int_t nColMax = digitsIn->GetNcol();
+ Int_t nRowMax = digitsIn->GetNrow();
Int_t nTimeTotal = digitsIn->GetNtime();
// Detector wise calibration object for the gain factors
- const AliTRDCalDet *calGainFactorDet = calibration->GetGainFactorDet();
+ const AliTRDCalDet *calGainFactorDet = calibration->GetGainFactorDet();
// Calibration object with pad wise values for the gain factors
- AliTRDCalROC *calGainFactorROC = calibration->GetGainFactorROC(idet);
+ AliTRDCalROC *calGainFactorROC = calibration->GetGainFactorROC(idet);
// Calibration value for chamber wise gain factor
- Float_t calGainFactorDetValue = calGainFactorDet->GetValue(idet);
+ Float_t calGainFactorDetValue = calGainFactorDet->GetValue(idet);
+
+
+ // Detector wise calibration object for the noise
+ const AliTRDCalDet *calNoiseDet = calibration->GetNoiseDet();
+ // Calibration object with pad wise values for the noise
+ AliTRDCalROC *calNoiseROC = calibration->GetNoiseROC(idet);
+ // Calibration value for chamber wise noise
+ Float_t calNoiseDetValue = calNoiseDet->GetValue(idet);
Int_t nClusters = 0;
- AliTRDdataArrayF *digitsOut = new AliTRDdataArrayF(digitsIn->GetNrow()
- ,digitsIn->GetNcol()
- ,digitsIn->GetNtime());
+ AliTRDdataArrayF *digitsOut = new AliTRDdataArrayF(nRowMax, nColMax, nTimeTotal);
+ AliTRDdataArrayS padStatus(nRowMax, nColMax, nTimeTotal);
ResetHelperIndexes(indexesIn);
// Apply the gain and the tail cancelation via digital filter
TailCancelation(digitsIn
- ,digitsOut
- ,indexesIn
- ,fIndexesOut
- ,nTimeTotal
+ ,digitsOut
+ ,indexesIn
+ ,fIndexesOut
+ ,nTimeTotal
,adcThreshold
,calGainFactorROC
,calGainFactorDetValue);
Int_t time = 0;
Int_t iPad = 0;
+ UChar_t status[3]={0, 0, 0}, ipos = 0;
fIndexesOut->ResetCounters();
- while (fIndexesOut->NextRCTbinIndex(row, col, time))
- {
+ while (fIndexesOut->NextRCTbinIndex(row, col, time)){
- Float_t signalM = TMath::Abs(digitsOut->GetDataUnchecked(row,col,time));
-
- // Look for the maximum
- if (signalM >= maxThresh)
- {
-
- if (col + 1 >= nColMax || col-1 < 0)
- continue;
+ Float_t signalM = TMath::Abs(digitsOut->GetDataUnchecked(row,col,time));
+ status[1] = digitsIn->GetPadStatus(row,col,time);
+ if(status[1]) SETBIT(ipos, AliTRDcluster::kMaskedCenter);
- Float_t signalL = TMath::Abs(digitsOut->GetDataUnchecked(row,col+1,time));
- Float_t signalR = TMath::Abs(digitsOut->GetDataUnchecked(row,col-1,time));
+ if(signalM < maxThresh) continue;
- if ((TMath::Abs(signalL) <= signalM) &&
- (TMath::Abs(signalR) < signalM))
- {
- if ((TMath::Abs(signalL) >= sigThresh) ||
- (TMath::Abs(signalR) >= sigThresh))
- {
- // Maximum found, mark the position by a negative signal
- digitsOut->SetDataUnchecked(row,col,time,-signalM);
- fIndexesMaxima->AddIndexTBin(row,col,time);
- }
- }
+ Float_t noiseMiddleThresh = minMaxCutSigma*calNoiseDetValue*calNoiseROC->GetValue(col,row);
+ if (signalM < noiseMiddleThresh) continue;
+ if (col + 1 >= nColMax || col-1 < 0) continue;
+
+ Float_t signalL = TMath::Abs(digitsOut->GetDataUnchecked(row,col+1,time));
+ status[0] = digitsIn->GetPadStatus(row,col+1,time);
+ if(status[0]) SETBIT(ipos, AliTRDcluster::kMaskedLeft);
+
+ Float_t signalR = TMath::Abs(digitsOut->GetDataUnchecked(row,col-1,time));
+ status[2] = digitsIn->GetPadStatus(row,col-1,time);
+ if(status[2]) SETBIT(ipos, AliTRDcluster::kMaskedRight);
+
+ // reject candidates with more than 1 problematic pad
+ if(ipos == 3 || ipos > 4) continue;
+
+ if(!status[1]){ // good central pad
+ if(!ipos){ // all pads are OK
+ if ((signalL <= signalM) && (signalR < signalM)) {
+ if ((signalL >= sigThresh) || (signalR >= sigThresh)) {
+ Float_t noiseSumThresh = minLeftRightCutSigma*calNoiseDetValue*calNoiseROC->GetValue(col,row);
+ if((signalL+signalR+signalM) >= noiseSumThresh){
+ // Maximum found, mark the position by a negative signal
+ digitsOut->SetDataUnchecked(row,col,time,-signalM);
+ fIndexesMaxima->AddIndexTBin(row,col,time);
+ padStatus.SetDataUnchecked(row, col, time, ipos);
+ }
+ }
}
-
+ } else { // one of the neighbouring pads are bad
+ if(status[0] && signalR < signalM && signalR >= sigThresh){
+ digitsOut->SetDataUnchecked(row,col,time,-signalM);
+ digitsOut->SetDataUnchecked(row, col, time+1, 0.);
+ fIndexesMaxima->AddIndexTBin(row,col,time);
+ padStatus.SetDataUnchecked(row, col, time, ipos);
+ } else if(status[2] && signalL <= signalM && signalL >= sigThresh){
+ digitsOut->SetDataUnchecked(row,col,time,-signalM);
+ digitsOut->SetDataUnchecked(row, col, time-1, 0.);
+ fIndexesMaxima->AddIndexTBin(row,col,time);
+ padStatus.SetDataUnchecked(row, col, time, ipos);
+ }
+ }
+ } else { // wrong maximum pad
+ if ((signalL >= sigThresh) || (signalR >= sigThresh)) {
+ // Maximum found, mark the position by a negative signal
+ digitsOut->SetDataUnchecked(row,col,time,-maxThresh);
+ fIndexesMaxima->AddIndexTBin(row,col,time);
+ padStatus.SetDataUnchecked(row, col, time, ipos);
+ }
}
-
- // The index to the first cluster of a given ROC
+ }
+
+ // The index to the first cluster of a given ROC
Int_t firstClusterROC = -1;
- // The number of cluster in a given ROC
- Int_t nClusterROC = 0;
+ // The number of cluster in a given ROC
+ Int_t nClusterROC = 0;
- // Now check the maxima and calculate the cluster position
- fIndexesMaxima->ResetCounters();
- while (fIndexesMaxima->NextRCTbinIndex(row, col, time))
- {
+ // Now check the maxima and calculate the cluster position
+ fIndexesMaxima->ResetCounters();
+ while (fIndexesMaxima->NextRCTbinIndex(row, col, time)) {
// Maximum found ?
- if (digitsOut->GetDataUnchecked(row,col,time) < 0.0)
- {
+ if (digitsOut->GetDataUnchecked(row,col,time) < 0.0) {
- for (iPad = 0; iPad < kNclus; iPad++)
- {
- Int_t iPadCol = col - 1 + iPad;
- clusterSignal[iPad] = TMath::Abs(digitsOut->GetDataUnchecked(row,iPadCol,time));
- }
+ for (iPad = 0; iPad < kNclus; iPad++) {
+ Int_t iPadCol = col - 1 + iPad;
+ clusterSignal[iPad] = TMath::Abs(digitsOut->GetDataUnchecked(row,iPadCol,time));
+ }
- // Count the number of pads in the cluster
- Int_t nPadCount = 0;
- Int_t ii;
- // Look to the left
- ii = 0;
- while (TMath::Abs(digitsOut->GetDataUnchecked(row,col-ii ,time)) >= sigThresh)
- {
- nPadCount++;
- ii++;
- if (col-ii < 0) break;
- }
- // Look to the right
- ii = 0;
- while (TMath::Abs(digitsOut->GetDataUnchecked(row,col+ii+1,time)) >= sigThresh)
- {
- nPadCount++;
- ii++;
- if (col+ii+1 >= nColMax) break;
- }
- nClusters++;
-
- // Look for 5 pad cluster with minimum in the middle
- Bool_t fivePadCluster = kFALSE;
- if (col < (nColMax - 3))
- {
- if (digitsOut->GetDataUnchecked(row,col+2,time) < 0)
- {
- fivePadCluster = kTRUE;
- }
- if ((fivePadCluster) && (col < (nColMax - 5)))
- {
- if (digitsOut->GetDataUnchecked(row,col+4,time) >= sigThresh)
- {
- fivePadCluster = kFALSE;
- }
- }
- if ((fivePadCluster) && (col > 1))
- {
- if (digitsOut->GetDataUnchecked(row,col-2,time) >= sigThresh)
- {
- fivePadCluster = kFALSE;
- }
- }
- }
+ // Count the number of pads in the cluster
+ Int_t nPadCount = 0;
+ Int_t ii;
+ // Look to the right
+ ii = 0;
+ while (TMath::Abs(digitsOut->GetDataUnchecked(row,col-ii ,time)) >= sigThresh) {
+ nPadCount++;
+ ii++;
+ if (col-ii < 0) break;
+ }
+ // Look to the left
+ ii = 0;
+ while (TMath::Abs(digitsOut->GetDataUnchecked(row,col+ii+1,time)) >= sigThresh){
+ nPadCount++;
+ ii++;
+ if (col+ii+1 >= nColMax) break;
+ }
+ nClusters++;
+
+ // Look for 5 pad cluster with minimum in the middle
+ Bool_t fivePadCluster = kFALSE;
+ if (col < (nColMax - 3)){
+ if (digitsOut->GetDataUnchecked(row,col+2,time) < 0) {
+ fivePadCluster = kTRUE;
+ }
+ if ((fivePadCluster) && (col < (nColMax - 5))) {
+ if (digitsOut->GetDataUnchecked(row,col+4,time) >= sigThresh){
+ fivePadCluster = kFALSE;
+ }
+ }
+ if ((fivePadCluster) && (col > 1)){
+ if (digitsOut->GetDataUnchecked(row,col-2,time) >= sigThresh){
+ fivePadCluster = kFALSE;
+ }
+ }
+ }
- // 5 pad cluster
- // Modify the signal of the overlapping pad for the left part
- // of the cluster which remains from a previous unfolding
- if (iUnfold)
- {
- clusterSignal[0] *= ratioLeft;
- iUnfold = 0;
- }
+ // 5 pad cluster
+ // Modify the signal of the overlapping pad for the left part
+ // of the cluster which remains from a previous unfolding
+ if (iUnfold) {
+ clusterSignal[0] *= ratioLeft;
+ iUnfold = 0;
+ }
- // Unfold the 5 pad cluster
- if (fivePadCluster)
- {
- for (iPad = 0; iPad < kNsig; iPad++)
- {
- padSignal[iPad] = TMath::Abs(digitsOut->GetDataUnchecked(row
- ,col-1+iPad
- ,time));
- }
- // Unfold the two maxima and set the signal on
- // the overlapping pad to the ratio
- ratioRight = Unfold(kEpsilon,iplan,padSignal);
- ratioLeft = 1.0 - ratioRight;
- clusterSignal[2] *= ratioRight;
- iUnfold = 1;
- }
+ // Unfold the 5 pad cluster
+ if (fivePadCluster){
+ for (iPad = 0; iPad < kNsig; iPad++) {
+ padSignal[iPad] = TMath::Abs(digitsOut->GetDataUnchecked(row
+ ,col-1+iPad
+ ,time));
+ }
+ // Unfold the two maxima and set the signal on
+ // the overlapping pad to the ratio
+ ratioRight = Unfold(kEpsilon,ilayer,padSignal);
+ ratioLeft = 1.0 - ratioRight;
+ clusterSignal[2] *= ratioRight;
+ iUnfold = 1;
+ }
- // The position of the cluster in COL direction relative to the center pad (pad units)
- Double_t clusterPosCol = 0.0;
- if (AliTRDReconstructor::RecoParam()->LUTOn())
- {
- // Calculate the position of the cluster by using the
- // lookup table method
- clusterPosCol = LUTposition(iplan,clusterSignal[0]
- ,clusterSignal[1]
- ,clusterSignal[2]);
- }
- else
- {
- // Calculate the position of the cluster by using the
- // center of gravity method
- for (Int_t i = 0; i < kNsig; i++)
- {
- padSignal[i] = 0.0;
- }
- padSignal[2] = TMath::Abs(digitsOut->GetDataUnchecked(row,col ,time)); // Central pad
- padSignal[1] = TMath::Abs(digitsOut->GetDataUnchecked(row,col-1,time)); // Left pad
- padSignal[3] = TMath::Abs(digitsOut->GetDataUnchecked(row,col+1,time)); // Right pad
- if ((col > 2) &&
- (TMath::Abs(digitsOut->GetDataUnchecked(row,col-2,time)) < padSignal[1]))
- {
- padSignal[0] = TMath::Abs(digitsOut->GetDataUnchecked(row,col-2,time));
- }
- if ((col < nColMax - 3) &&
- (TMath::Abs(digitsOut->GetDataUnchecked(row,col+2,time)) < padSignal[3]))
- {
- padSignal[4] = TMath::Abs(digitsOut->GetDataUnchecked(row,col+2,time));
- }
- clusterPosCol = GetCOG(padSignal);
- }
+ // The position of the cluster in COL direction relative to the center pad (pad units)
+ Double_t clusterPosCol = 0.0;
+ if (fReconstructor->GetRecoParam() ->IsLUT()) {
+ // Calculate the position of the cluster by using the
+ // lookup table method
+ clusterPosCol = LUTposition(ilayer,clusterSignal[2]
+ ,clusterSignal[1]
+ ,clusterSignal[0]);
+ }
+ else {
+ // Calculate the position of the cluster by using the
+ // center of gravity method
+ for (Int_t i = 0; i < kNsig; i++) {
+ padSignal[i] = 0.0;
+ }
+ padSignal[2] = TMath::Abs(digitsOut->GetDataUnchecked(row,col ,time)); // Central pad
+ padSignal[1] = TMath::Abs(digitsOut->GetDataUnchecked(row,col+1,time)); // Left pad
+ padSignal[3] = TMath::Abs(digitsOut->GetDataUnchecked(row,col-1,time)); // Right pad
+ if ((col > 2) &&
+ (TMath::Abs(digitsOut->GetDataUnchecked(row,col-2,time)) < padSignal[1])) {
+ padSignal[4] = TMath::Abs(digitsOut->GetDataUnchecked(row,col-2,time));
+ }
+ if ((col < nColMax - 3) &&
+ (TMath::Abs(digitsOut->GetDataUnchecked(row,col+2,time)) < padSignal[3])){
+ padSignal[0] = TMath::Abs(digitsOut->GetDataUnchecked(row,col+2,time));
+ }
+ clusterPosCol = GetCOG(padSignal);
+ }
- // Store the amplitudes of the pads in the cluster for later analysis
- Short_t signals[7] = { 0, 0, 0, 0, 0, 0, 0 };
- for (Int_t jPad = col-3; jPad <= col+3; jPad++)
- {
- if ((jPad < 0) ||
- (jPad >= nColMax-1))
- {
- continue;
- }
- signals[jPad-col+3] = TMath::Nint(TMath::Abs(digitsOut->GetDataUnchecked(row,jPad,time)));
- }
+ // Store the amplitudes of the pads in the cluster for later analysis
+ // and check whether one of these pads is masked in the database
+ Short_t signals[7] = { 0, 0, 0, 0, 0, 0, 0 };
+ for (Int_t jPad = col-3; jPad <= col+3; jPad++) {
+ if ((jPad < 0) ||
+ (jPad >= nColMax-1)) {
+ continue;
+ }
+ signals[jPad-col+3] = TMath::Nint(TMath::Abs(digitsOut->GetDataUnchecked(row,jPad,time)));
+ }
- // Transform the local cluster coordinates into calibrated
- // space point positions defined in the local tracking system.
- // Here the calibration for T0, Vdrift and ExB is applied as well.
- Double_t clusterXYZ[6];
- clusterXYZ[0] = clusterPosCol;
- clusterXYZ[1] = clusterSignal[0];
- clusterXYZ[2] = clusterSignal[1];
- clusterXYZ[3] = clusterSignal[2];
- clusterXYZ[4] = 0.0;
- clusterXYZ[5] = 0.0;
- Int_t clusterRCT[3];
- clusterRCT[0] = row;
- clusterRCT[1] = col;
- clusterRCT[2] = 0;
-
- Bool_t out = kTRUE;
- if (fTransform->Transform(clusterXYZ, clusterRCT, ((UInt_t) time), out, 0)) {
-
- // Add the cluster to the output array
- // The track indices will be stored later
- Float_t clusterPos[3];
- clusterPos[0] = clusterXYZ[0];
- clusterPos[1] = clusterXYZ[1];
- clusterPos[2] = clusterXYZ[2];
- Float_t clusterSig[2];
- clusterSig[0] = clusterXYZ[4];
- clusterSig[1] = clusterXYZ[5];
- Double_t clusterCharge = clusterXYZ[3];
- Char_t clusterTimeBin = ((Char_t) clusterRCT[2]);
- AliTRDcluster *cluster = new AliTRDcluster(idet
- ,clusterCharge
- ,clusterPos
- ,clusterSig
- ,0x0
- ,((Char_t) nPadCount)
- ,signals
- ,((UChar_t) col)
- ,((UChar_t) row)
- ,((UChar_t) time)
- ,clusterTimeBin
- ,clusterPosCol
- ,volid);
- cluster->SetInChamber(!out);
-
- // Temporarily store the row, column and time bin of the center pad
- // Used to later on assign the track indices
- cluster->SetLabel( row,0);
- cluster->SetLabel( col,1);
- cluster->SetLabel(time,2);
-
- RecPoints()->Add(cluster);
-
- // Store the index of the first cluster in the current ROC
- if (firstClusterROC < 0)
- {
- firstClusterROC = RecPoints()->GetEntriesFast() - 1;
- }
-
- // Count the number of cluster in the current ROC
- nClusterROC++;
-
- } // if: Transform ok ?
-
- } // if: Maximum found ?
+ // Transform the local cluster coordinates into calibrated
+ // space point positions defined in the local tracking system.
+ // Here the calibration for T0, Vdrift and ExB is applied as well.
+ Double_t clusterXYZ[6];
+ clusterXYZ[0] = clusterPosCol;
+ clusterXYZ[1] = clusterSignal[2];
+ clusterXYZ[2] = clusterSignal[1];
+ clusterXYZ[3] = clusterSignal[0];
+ clusterXYZ[4] = 0.0;
+ clusterXYZ[5] = 0.0;
+ Int_t clusterRCT[3];
+ clusterRCT[0] = row;
+ clusterRCT[1] = col;
+ clusterRCT[2] = 0;
+
+ Bool_t out = kTRUE;
+ if (fTransform->Transform(clusterXYZ, clusterRCT, ((UInt_t) time), out, 0)) {
+
+ // Add the cluster to the output array
+ // The track indices will be stored later
+ Float_t clusterPos[3];
+ clusterPos[0] = clusterXYZ[0];
+ clusterPos[1] = clusterXYZ[1];
+ clusterPos[2] = clusterXYZ[2];
+ Float_t clusterSig[2];
+ clusterSig[0] = clusterXYZ[4];
+ clusterSig[1] = clusterXYZ[5];
+ Double_t clusterCharge = clusterXYZ[3];
+ Char_t clusterTimeBin = ((Char_t) clusterRCT[2]);
+ AliTRDcluster *cluster = new AliTRDcluster(idet
+ ,clusterCharge
+ ,clusterPos
+ ,clusterSig
+ ,0x0
+ ,((Char_t) nPadCount)
+ ,signals
+ ,((UChar_t) col)
+ ,((UChar_t) row)
+ ,((UChar_t) time)
+ ,clusterTimeBin
+ ,clusterPosCol
+ ,volid);
+ cluster->SetInChamber(!out);
+
+ UChar_t maskPosition = padStatus.GetDataUnchecked(row, col, time);
+ if (maskPosition) {
+ cluster->SetPadMaskedPosition(maskPosition);
+ if (maskPosition & AliTRDcluster::kMaskedLeft) {
+ cluster->SetPadMaskedStatus(status[0]);
+ }
+ else if (maskPosition & AliTRDcluster::kMaskedCenter) {
+ cluster->SetPadMaskedStatus(status[1]);
+ }
+ else {
+ cluster->SetPadMaskedStatus(status[2]);
+ }
+ }
- }
+ // Temporarily store the row, column and time bin of the center pad
+ // Used to later on assign the track indices
+ cluster->SetLabel( row,0);
+ cluster->SetLabel( col,1);
+ cluster->SetLabel(time,2);
+
+ RecPoints()->Add(cluster);
+
+ // Store the index of the first cluster in the current ROC
+ if (firstClusterROC < 0) {
+ firstClusterROC = RecPoints()->GetEntriesFast() - 1;
+ }
+
+ // Count the number of cluster in the current ROC
+ nClusterROC++;
+
+ } // if: Transform ok ?
+ } // if: Maximum found ?
+ }
delete digitsOut;
- if (fAddLabels)
- {
- AddLabels(idet, firstClusterROC, nClusterROC);
- }
+ if (fAddLabels) AddLabels(idet, firstClusterROC, nClusterROC);
// Write the cluster and reset the array
WriteClusters(idet);
+ signal[3]
+ signal[4];
+ // ???????????? CBL
Double_t res = (0.0 * (-signal[0] + signal[4])
+ (-signal[1] + signal[3])) / sum;
}
//_____________________________________________________________________________
-Double_t AliTRDclusterizer::Unfold(Double_t eps, Int_t plane, Double_t *padSignal)
+Double_t AliTRDclusterizer::Unfold(Double_t eps, Int_t layer, Double_t *padSignal)
{
//
// Method to unfold neighbouring maxima.
/ ((1.0 - ratio)*padSignal[2] + padSignal[3] + padSignal[4]);
// Set cluster charge ratio
- irc = calibration->PadResponse(1.0,maxLeft ,plane,newSignal);
+ irc = calibration->PadResponse(1.0,maxLeft ,layer,newSignal);
Double_t ampLeft = padSignal[1] / newSignal[1];
- irc = calibration->PadResponse(1.0,maxRight,plane,newSignal);
+ irc = calibration->PadResponse(1.0,maxRight,layer,newSignal);
Double_t ampRight = padSignal[3] / newSignal[1];
// Apply pad response to parameters
- irc = calibration->PadResponse(ampLeft ,maxLeft ,plane,newLeftSignal );
- irc = calibration->PadResponse(ampRight,maxRight,plane,newRightSignal);
+ irc = calibration->PadResponse(ampLeft ,maxLeft ,layer,newLeftSignal );
+ irc = calibration->PadResponse(ampRight,maxRight,layer,newRightSignal);
// Calculate new overlapping ratio
ratio = TMath::Min((Double_t) 1.0
}
//_____________________________________________________________________________
-void AliTRDclusterizer::TailCancelation(AliTRDdataArrayS *digitsIn
+void AliTRDclusterizer::TailCancelation(AliTRDdataArrayDigits *digitsIn
, AliTRDdataArrayF *digitsOut
, AliTRDSignalIndex *indexesIn
, AliTRDSignalIndex *indexesOut
Double_t gain = calGainFactorDetValue
* calGainFactorROCValue;
+ Bool_t corrupted = kFALSE;
for (iTime = 0; iTime < nTimeTotal; iTime++)
{
// Apply gain gain factor
inADC[iTime] = digitsIn->GetDataUnchecked(iRow,iCol,iTime);
+ if(digitsIn->GetPadStatus(iRow, iCol, iTime)) corrupted = kTRUE;
inADC[iTime] /= gain;
outADC[iTime] = inADC[iTime];
}
-
- // Apply the tail cancelation via the digital filter
- if (AliTRDReconstructor::RecoParam()->TCOn())
- {
- DeConvExp(inADC,outADC,nTimeTotal,AliTRDReconstructor::RecoParam()->GetTCnexp());
- }
+ if(!corrupted)
+ {
+ // Apply the tail cancelation via the digital filter
+ // (only for non-coorupted pads)
+ if (fReconstructor->GetRecoParam() ->IsTailCancelation())
+ {
+ DeConvExp(inADC,outADC,nTimeTotal,fReconstructor->GetRecoParam() ->GetTCnexp());
+ }
+ }
indexesIn->ResetTbinCounter();
while (indexesIn->NextTbinIndex(iTime))
c2 = 0.000;
}
if (nexp == 2) { // 2 Exponentials
- r1 = 1.156;
- r2 = 0.130;
- c1 = 0.114;
- c2 = 0.624;
+ Double_t par[4];
+ fReconstructor->GetTCParams(par);
+ r1 = par[0];//1.156;
+ r2 = par[1];//0.130;
+ c1 = par[2];//0.114;
+ c2 = par[3];//0.624;
}
coefficients[0] = c1;
const Int_t kNlut = 128;
- fLUTbin = AliTRDgeometry::kNplan * kNlut;
+ fLUTbin = AliTRDgeometry::kNlayer * kNlut;
// The lookup table from Bogdan
- Float_t lut[AliTRDgeometry::kNplan][kNlut] = {
+ Float_t lut[AliTRDgeometry::kNlayer][kNlut] = {
{
0.0070, 0.0150, 0.0224, 0.0298, 0.0374, 0.0454, 0.0533, 0.0611,
0.0684, 0.0755, 0.0827, 0.0900, 0.0975, 0.1049, 0.1120, 0.1187,
}
fLUT = new Double_t[fLUTbin];
- for (Int_t iplan = 0; iplan < AliTRDgeometry::kNplan; iplan++) {
+ for (Int_t ilayer = 0; ilayer < AliTRDgeometry::kNlayer; ilayer++) {
for (Int_t ilut = 0; ilut < kNlut; ilut++ ) {
- fLUT[iplan*kNlut+ilut] = lut[iplan][ilut];
+ fLUT[ilayer*kNlut+ilut] = lut[ilayer][ilut];
}
}
}
//_____________________________________________________________________________
-Double_t AliTRDclusterizer::LUTposition(Int_t iplane, Double_t ampL
+Double_t AliTRDclusterizer::LUTposition(Int_t ilayer, Double_t ampL
, Double_t ampC, Double_t ampR) const
{
//
Int_t side = 0;
Int_t ix;
- Double_t xMin[AliTRDgeometry::kNplan] = { 0.006492, 0.006377, 0.006258
- , 0.006144, 0.006030, 0.005980 };
- Double_t xMax[AliTRDgeometry::kNplan] = { 0.960351, 0.965870, 0.970445
- , 0.974352, 0.977667, 0.996101 };
+ Double_t xMin[AliTRDgeometry::kNlayer] = { 0.006492, 0.006377, 0.006258
+ , 0.006144, 0.006030, 0.005980 };
+ Double_t xMax[AliTRDgeometry::kNlayer] = { 0.960351, 0.965870, 0.970445
+ , 0.974352, 0.977667, 0.996101 };
if (ampL > ampR) {
x = (ampL - ampR) / ampC;
if (ampL != ampR) {
- xmin = xMin[iplane] + 0.000005;
- xmax = xMax[iplane] - 0.000005;
+ xmin = xMin[ilayer] + 0.000005;
+ xmax = xMax[ilayer] - 0.000005;
xwid = (xmax - xmin) / 127.0;
if (x < xmin) {
}
else {
ix = (Int_t) ((x - xmin) / xwid);
- pos = side * fLUT[iplane*kNlut+ix];
+ pos = side * fLUT[ilayer*kNlut+ix];
}
}