[u/mrichter/AliRoot.git] / HLT / src / AliL3ClustFinderNew.cxx

//$Id$

// Author: Anders Vestbo <mailto:vestbo@fi.uib.no>
//*-- Copyright &copy ASV 


#include "AliL3StandardIncludes.h"

#include "AliL3Logging.h"
#include "AliL3ClustFinderNew.h"
#include "AliL3DigitData.h"
#include "AliL3Transform.h"
#include "AliL3SpacePointData.h"
#include "AliL3MemHandler.h"

#if GCCVERSION == 3
using namespace std;
#endif

/** \class AliL3ClustFinderNew
//<pre>
//_____________________________________________________________
// AliL3ClustFinderNew
//
// The current cluster finder for HLT
// Based on STAR L3
//</pre> */

ClassImp(AliL3ClustFinderNew)

AliL3ClustFinderNew::AliL3ClustFinderNew()
{
  fMatch = 2;
  fThreshold = 10;
  fXYErr = 0.2;
  fZErr = 0.3;
  fDeconvPad = kTRUE;
  fDeconvTime = kTRUE;
  fstdout = kFALSE;
  fcalcerr = kTRUE;
}

AliL3ClustFinderNew::~AliL3ClustFinderNew()
{
}

void AliL3ClustFinderNew::InitSlice(Int_t slice,Int_t patch,Int_t firstrow, Int_t lastrow,Int_t nmaxpoints)
{
  fNClusters = 0;
  fMaxNClusters = nmaxpoints;
  fCurrentSlice = slice;
  fCurrentPatch = patch;
  fFirstRow = firstrow;
  fLastRow = lastrow;
}

void AliL3ClustFinderNew::InitSlice(Int_t slice,Int_t patch,Int_t nmaxpoints)
{
  fNClusters = 0;
  fMaxNClusters = nmaxpoints;
  fCurrentSlice = slice;
  fCurrentPatch = patch;
}

void AliL3ClustFinderNew::SetOutputArray(AliL3SpacePointData *pt)
{
  fSpacePointData = pt;
}

void AliL3ClustFinderNew::Read(UInt_t ndigits,AliL3DigitRowData *ptr)
{
  fNDigitRowData = ndigits;
  fDigitRowData = ptr;
}

void AliL3ClustFinderNew::ProcessDigits()
{
  //Loop over rows, and call processrow
  
  
  AliL3DigitRowData *tempPt = (AliL3DigitRowData*)fDigitRowData;
  
  for(Int_t i=fFirstRow; i<=fLastRow; i++)
    {
      fCurrentRow = i;
      ProcessRow(tempPt);
      Byte_t *tmp = (Byte_t*)tempPt;
      Int_t size = sizeof(AliL3DigitRowData) + tempPt->fNDigit*sizeof(AliL3DigitData);
      tmp += size;
      tempPt = (AliL3DigitRowData*)tmp;
    }
  LOG(AliL3Log::kInformational,"AliL3ClustFinderNew::WriteClusters","Space points")
    <<"Cluster finder found "<<fNClusters<<" clusters in slice "<<fCurrentSlice
    <<" patch "<<fCurrentPatch<<ENDLOG;
}

void AliL3ClustFinderNew::ProcessRow(AliL3DigitRowData *tempPt)
{

  UInt_t last_pad = 123456789;

  ClusterData *pad1[2500]; //2 lists for internal memory=2pads
  ClusterData *pad2[2500]; //2 lists for internal memory=2pads
  ClusterData clusterlist[5000]; //Clusterlist

  ClusterData **currentPt;  //List of pointers to the current pad
  ClusterData **previousPt; //List of pointers to the previous pad
  currentPt = pad2;
  previousPt = pad1;
  UInt_t n_previous=0,n_current=0,n_total=0;

  //Loop over sequences of this row:
  for(UInt_t bin=0; bin<tempPt->fNDigit; bin++)
    {
      AliL3DigitData *data = tempPt->fDigitData;
      if(data[bin].fPad != last_pad)
	{
	  //This is a new pad
	  
	  //Switch:
	  if(currentPt == pad2)
	    {
	      currentPt = pad1;
	      previousPt = pad2;
	    }
	  else 
	    {
	      currentPt = pad2;
	      previousPt = pad1;
	    }
	  n_previous = n_current;
	  n_current = 0;
	  if(bin[data].fPad != last_pad+1)
	    {
	      //this happens if there is a pad with no signal.
	      n_previous = n_current = 0;
	    }
	  last_pad = data[bin].fPad;
	}

      Bool_t new_cluster = kTRUE;
      UInt_t seq_charge=0,seq_average=0,seq_error=0;
      UInt_t last_charge=0,last_was_falling=0;
      Int_t new_bin=-1;

      if(fDeconvTime)
	{
	redo: //This is a goto.
	  if(new_bin > -1)
	    {
	      bin = new_bin;
	      new_bin = -1;
	    }
	  
	  last_charge=0;
	  last_was_falling = 0;
	}
      
      while(1) //Loop over current sequence
	{
	  if(data[bin].fTime >= AliL3Transform::GetNTimeBins())
	    {
	      LOG(AliL3Log::kFatal,"AliL3ClustFinderNew::ProcessRow","Digits")
		<<"Timebin out of range "<<(Int_t)data[bin].fTime<<ENDLOG;
	      break;
	    }

	  //Get the current ADC-value
	  UInt_t charge = data[bin].fCharge;
	  
	  if(fDeconvTime)
	    {
	      //Check if the last pixel in the sequence is smaller than this
	      if(charge > last_charge)
		{
		  if(last_was_falling)
		    {
		      new_bin = bin;
		      break;
		    }
		}
	      else last_was_falling = 1; //last pixel was larger than this
	      last_charge = charge;
	    }
	  
	  //Sum the total charge of this sequence
	  seq_charge += charge;
	  seq_average += data[bin].fTime*charge;
	  seq_error += data[bin].fTime*data[bin].fTime*charge;

	  //Check where to stop:
	  if(data[bin+1].fPad != data[bin].fPad) //new pad
	    break; 
	  if(data[bin+1].fTime != data[bin].fTime+1) //end of sequence
	    break;
	  
	  bin++;
	}//end loop over sequence
      
      //Calculate mean of sequence:
      Int_t seq_mean=0;
      if(seq_charge)
	seq_mean = seq_average/seq_charge;
      else
	{
	  LOG(AliL3Log::kFatal,"AliL3ClustFinderNew::ProcessRow","Data")
	    <<"Error in data given to the cluster finder"<<ENDLOG;
	  seq_mean = 1;
	  seq_charge = 1;
	}
      
      //Calculate mean in pad direction:
      Int_t pad_mean = seq_charge*data[bin].fPad;
      Int_t pad_error = data[bin].fPad*pad_mean;
      
      //Compare with results on previous pad:
      for(UInt_t p=0; p<n_previous; p++)
	{
	  //dont merge sequences on the same pad twice
	  if(previousPt[p]->fLastMergedPad==data[bin].fPad) continue;

	  Int_t difference = seq_mean - previousPt[p]->fMean;
	  if(difference < -fMatch) break;

	  if(difference <= fMatch) //There is a match here!!
	    {
	      ClusterData *local = previousPt[p];

	      if(fDeconvPad)
		{
		  if(seq_charge > local->fLastCharge)
		    {
		      if(local->fChargeFalling) //The previous pad was falling
			{			
			  break; //create a new cluster
			}		    
		    }
		  else
		    local->fChargeFalling = 1;
		  local->fLastCharge = seq_charge;
		}
	      
	      //Don't create a new cluster, because we found a match
	      new_cluster = kFALSE;
	      
	      //Update cluster on current pad with the matching one:
	      local->fTotalCharge += seq_charge;
	      local->fPad += pad_mean;
	      local->fPad2 += pad_error;
	      local->fTime += seq_average;
	      local->fTime2 += seq_error;
	      local->fMean = seq_mean;
	      local->fFlags++; //means we have more than one pad 
	      local->fLastMergedPad = data[bin].fPad;

	      currentPt[n_current] = local;
	      n_current++;
	      
	      break;
	    } //Checking for match at previous pad
	} //Loop over results on previous pad.
      
      if(new_cluster)
	{
	  //Start a new cluster. Add it to the clusterlist, and update
	  //the list of pointers to clusters in current pad.
	  //current pad will be previous pad on next pad.

	  //Add to the clusterlist:
	  ClusterData *tmp = &clusterlist[n_total];
	  tmp->fTotalCharge = seq_charge;
	  tmp->fPad = pad_mean;
	  tmp->fPad2 = pad_error;
	  tmp->fTime = seq_average;
	  tmp->fTime2 = seq_error;
	  tmp->fMean = seq_mean;
	  tmp->fFlags = 0;  //flags for 1 pad clusters
	  tmp->fLastMergedPad = data[bin].fPad;
	  if(fDeconvPad)
	    {
	      tmp->fChargeFalling = 0;
	      tmp->fLastCharge = seq_charge;
	    }

	  //Update list of pointers to previous pad:
	  currentPt[n_current] = &clusterlist[n_total];
	  n_total++;
	  n_current++;
	}

      if(fDeconvTime)
	if(new_bin >= 0) goto redo;
    }//Loop over digits on this padrow
  
  WriteClusters(n_total,clusterlist);
}

void AliL3ClustFinderNew::WriteClusters(Int_t n_clusters,ClusterData *list)
{
  Int_t thisrow,thissector;
  UInt_t counter = fNClusters;
  
  for(int j=0; j<n_clusters; j++)
    {
      if(!list[j].fFlags) continue; //discard 1 pad clusters
      if(list[j].fTotalCharge < fThreshold) continue; //noise cluster

      Float_t xyz[3];      
      Float_t fpad =(Float_t)list[j].fPad /(Float_t)list[j].fTotalCharge;
      Float_t fpad2=fXYErr;
      Float_t ftime =(Float_t)list[j].fTime /(Float_t)list[j].fTotalCharge;
      Float_t ftime2=fZErr;

      if(fcalcerr) {
	fpad2=(Float_t)list[j].fPad2/(Float_t)list[j].fTotalCharge - fpad*fpad;
	fpad2 = sqrt(fpad2);
	ftime2=(Float_t)list[j].fTime2/(Float_t)list[j].fTotalCharge - ftime*ftime;
	ftime2 = sqrt(ftime2); 
      }
       
      if(fstdout==kTRUE)
	cout<<"WriteCluster: padrow "<<fCurrentRow<<" pad "<<fpad << " +- "<<fpad2<<" time "<<ftime<<" +- "<<ftime2<<" charge "<<list[j].fTotalCharge<<endl;

      AliL3Transform::Slice2Sector(fCurrentSlice,fCurrentRow,thissector,thisrow);
      AliL3Transform::Raw2Local(xyz,thissector,thisrow,fpad,ftime);
      if(xyz[0]==0) LOG(AliL3Log::kError,"AliL3ClustFinder","Cluster Finder")
	<<AliL3Log::kDec<<"Zero cluster"<<ENDLOG;
      if(fNClusters >= fMaxNClusters)
	{
	  LOG(AliL3Log::kError,"AliL3ClustFinder::WriteClusters","Cluster Finder")
	    <<AliL3Log::kDec<<"Too many clusters"<<ENDLOG;
	  return;
	}  
      fSpacePointData[counter].fCharge = list[j].fTotalCharge;
      fSpacePointData[counter].fX = xyz[0];
      fSpacePointData[counter].fY = xyz[1];
      fSpacePointData[counter].fZ = xyz[2];
      fSpacePointData[counter].fPadRow = fCurrentRow;
      fSpacePointData[counter].fXYErr = fpad2;
      fSpacePointData[counter].fZErr  = ftime2;
      fSpacePointData[counter].fID = counter
	+((fCurrentSlice&0x7f)<<25)+((fCurrentPatch&0x7)<<22);//Uli
#ifdef do_mc
      Int_t trackID[3];
      GetTrackID((Int_t)rint(fpad),(Int_t)rint(ftime),trackID);
      //cout<<"padrow "<<fCurrentRow<<" pad "<<(Int_t)rint(fpad)<<" time "<<(Int_t)rint(ftime)<<" Trackid "<<trackID[0]<<endl;
      fSpacePointData[counter].fTrackID[0] = trackID[0];
      fSpacePointData[counter].fTrackID[1] = trackID[1];
      fSpacePointData[counter].fTrackID[2] = trackID[2];
#endif

      fNClusters++;
      counter++;
    }
}

#ifdef do_mc
void AliL3ClustFinderNew::GetTrackID(Int_t pad,Int_t time,Int_t *trackID)
{
  AliL3DigitRowData *rowPt = (AliL3DigitRowData*)fDigitRowData;
  
  trackID[0]=trackID[1]=trackID[2]=-2;
  //cout<<"Looking for pad "<<pad<<" time "<<time<<endl;
  for(Int_t i=fFirstRow; i<=fLastRow; i++)
    {
      if(rowPt->fRow < (UInt_t)fCurrentRow)
	{
	  AliL3MemHandler::UpdateRowPointer(rowPt);
	  continue;
	}
      AliL3DigitData *digPt = (AliL3DigitData*)rowPt->fDigitData;
      for(UInt_t j=0; j<rowPt->fNDigit; j++)
	{
	  Int_t cpad = digPt[j].fPad;
	  Int_t ctime = digPt[j].fTime;
	  if(cpad != pad) continue;
	  if(ctime != time) continue;
	  //if(cpad != pad && ctime != ctime) continue;
	  //cout<<"Reading row "<<fCurrentRow<<" pad "<<cpad<<" time "<<ctime<<" trackID "<<digPt[j].fTrackID[0]<<endl;
	  trackID[0] = digPt[j].fTrackID[0];
	  trackID[1] = digPt[j].fTrackID[1];
	  trackID[2] = digPt[j].fTrackID[2];
	  break;
	  //cout<<"Reading trackID "<<trackID[0]<<endl;
	}
      break;
    }
  
}
#endif
Commit	Line	Data
4ec96209	1	//$Id$
4ec96209	2
b661165c	3	// Author: Anders Vestbo <mailto:vestbo@fi.uib.no>
b661165c	4	//*-- Copyright &copy ASV
46fbeb8a	5
118c26c3	6
	7	#include "AliL3StandardIncludes.h"
	8
46fbeb8a	9	#include "AliL3Logging.h"
	10	#include "AliL3ClustFinderNew.h"
	11	#include "AliL3DigitData.h"
	12	#include "AliL3Transform.h"
	13	#include "AliL3SpacePointData.h"
10f815d9	14	#include "AliL3MemHandler.h"
46fbeb8a	15
118c26c3	16	#if GCCVERSION == 3
	17	using namespace std;
	18	#endif
	19
b62fa469	20	/** \class AliL3ClustFinderNew
b62fa469	21	//<pre>
b661165c	22	//_____________________________________________________________
	23	// AliL3ClustFinderNew
	24	//
	25	// The current cluster finder for HLT
	26	// Based on STAR L3
b62fa469	27	//</pre> */
b661165c	28
46fbeb8a	29	ClassImp(AliL3ClustFinderNew)
	30
	31	AliL3ClustFinderNew::AliL3ClustFinderNew()
	32	{
9953e05f	33	fMatch = 2;
77f08cd3	34	fThreshold = 10;
	35	fXYErr = 0.2;
	36	fZErr = 0.3;
99a6d5c9	37	fDeconvPad = kTRUE;
2bbf57d8	38	fDeconvTime = kTRUE;
77f08cd3	39	fstdout = kFALSE;
c3dd27a3	40	fcalcerr = kTRUE;
46fbeb8a	41	}
46fbeb8a	42
46fbeb8a	43	AliL3ClustFinderNew::~AliL3ClustFinderNew()
46fbeb8a	44	{
46fbeb8a	45	}
	46
	47	void AliL3ClustFinderNew::InitSlice(Int_t slice,Int_t patch,Int_t firstrow, Int_t lastrow,Int_t nmaxpoints)
	48	{
	49	fNClusters = 0;
	50	fMaxNClusters = nmaxpoints;
	51	fCurrentSlice = slice;
	52	fCurrentPatch = patch;
	53	fFirstRow = firstrow;
	54	fLastRow = lastrow;
	55	}
	56
	57	void AliL3ClustFinderNew::InitSlice(Int_t slice,Int_t patch,Int_t nmaxpoints)
	58	{
	59	fNClusters = 0;
	60	fMaxNClusters = nmaxpoints;
	61	fCurrentSlice = slice;
	62	fCurrentPatch = patch;
	63	}
	64
	65	void AliL3ClustFinderNew::SetOutputArray(AliL3SpacePointData *pt)
	66	{
46fbeb8a	67	fSpacePointData = pt;
	68	}
	69
46fbeb8a	70	void AliL3ClustFinderNew::Read(UInt_t ndigits,AliL3DigitRowData *ptr)
	71	{
	72	fNDigitRowData = ndigits;
	73	fDigitRowData = ptr;
46fbeb8a	74	}
	75
	76	void AliL3ClustFinderNew::ProcessDigits()
	77	{
	78	//Loop over rows, and call processrow
	79
	80
	81	AliL3DigitRowData tempPt = (AliL3DigitRowData)fDigitRowData;
	82
	83	for(Int_t i=fFirstRow; i<=fLastRow; i++)
	84	{
	85	fCurrentRow = i;
	86	ProcessRow(tempPt);
	87	Byte_t tmp = (Byte_t)tempPt;
	88	Int_t size = sizeof(AliL3DigitRowData) + tempPt->fNDigit*sizeof(AliL3DigitData);
	89	tmp += size;
	90	tempPt = (AliL3DigitRowData*)tmp;
	91	}
355debe1	92	LOG(AliL3Log::kInformational,"AliL3ClustFinderNew::WriteClusters","Space points")
	93	<<"Cluster finder found "<<fNClusters<<" clusters in slice "<<fCurrentSlice
	94	<<" patch "<<fCurrentPatch<<ENDLOG;
46fbeb8a	95	}
46fbeb8a	96
46fbeb8a	97	void AliL3ClustFinderNew::ProcessRow(AliL3DigitRowData *tempPt)
	98	{
	99
	100	UInt_t last_pad = 123456789;
	101
	102	ClusterData *pad1[2500]; //2 lists for internal memory=2pads
	103	ClusterData *pad2[2500]; //2 lists for internal memory=2pads
	104	ClusterData clusterlist[5000]; //Clusterlist
	105
77f08cd3	106	ClusterData **currentPt; //List of pointers to the current pad
46fbeb8a	107	ClusterData **previousPt; //List of pointers to the previous pad
	108	currentPt = pad2;
	109	previousPt = pad1;
	110	UInt_t n_previous=0,n_current=0,n_total=0;
	111
	112	//Loop over sequences of this row:
	113	for(UInt_t bin=0; bin<tempPt->fNDigit; bin++)
	114	{
46fbeb8a	115	AliL3DigitData *data = tempPt->fDigitData;
	116	if(data[bin].fPad != last_pad)
	117	{
	118	//This is a new pad
46fbeb8a	119
	120	//Switch:
	121	if(currentPt == pad2)
	122	{
	123	currentPt = pad1;
	124	previousPt = pad2;
46fbeb8a	125	}
	126	else
	127	{
	128	currentPt = pad2;
	129	previousPt = pad1;
	130	}
	131	n_previous = n_current;
	132	n_current = 0;
99a6d5c9	133	if(bin[data].fPad != last_pad+1)
	134	{
	135	//this happens if there is a pad with no signal.
	136	n_previous = n_current = 0;
	137	}
	138	last_pad = data[bin].fPad;
46fbeb8a	139	}
	140
	141	Bool_t new_cluster = kTRUE;
c3dd27a3	142	UInt_t seq_charge=0,seq_average=0,seq_error=0;
99a6d5c9	143	UInt_t last_charge=0,last_was_falling=0;
99a6d5c9	144	Int_t new_bin=-1;
77f08cd3	145
99a6d5c9	146	if(fDeconvTime)
	147	{
	148	redo: //This is a goto.
	149	if(new_bin > -1)
	150	{
	151	bin = new_bin;
	152	new_bin = -1;
	153	}
	154
	155	last_charge=0;
	156	last_was_falling = 0;
	157	}
	158
46fbeb8a	159	while(1) //Loop over current sequence
46fbeb8a	160	{
494fad94	161	if(data[bin].fTime >= AliL3Transform::GetNTimeBins())
adb0e845	162	{
	163	LOG(AliL3Log::kFatal,"AliL3ClustFinderNew::ProcessRow","Digits")
	164	<<"Timebin out of range "<<(Int_t)data[bin].fTime<<ENDLOG;
	165	break;
	166	}
	167
46fbeb8a	168	//Get the current ADC-value
46fbeb8a	169	UInt_t charge = data[bin].fCharge;
99a6d5c9	170
	171	if(fDeconvTime)
	172	{
	173	//Check if the last pixel in the sequence is smaller than this
	174	if(charge > last_charge)
	175	{
	176	if(last_was_falling)
	177	{
	178	new_bin = bin;
	179	break;
	180	}
	181	}
	182	else last_was_falling = 1; //last pixel was larger than this
	183	last_charge = charge;
	184	}
	185
46fbeb8a	186	//Sum the total charge of this sequence
46fbeb8a	187	seq_charge += charge;
46fbeb8a	188	seq_average += data[bin].fTime*charge;
c3dd27a3	189	seq_error += data[bin].fTimedata[bin].fTimecharge;
c3dd27a3	190
46fbeb8a	191	//Check where to stop:
	192	if(data[bin+1].fPad != data[bin].fPad) //new pad
	193	break;
	194	if(data[bin+1].fTime != data[bin].fTime+1) //end of sequence
	195	break;
	196
	197	bin++;
	198	}//end loop over sequence
	199
	200	//Calculate mean of sequence:
	201	Int_t seq_mean=0;
	202	if(seq_charge)
	203	seq_mean = seq_average/seq_charge;
	204	else
99a6d5c9	205	{
adb0e845	206	LOG(AliL3Log::kFatal,"AliL3ClustFinderNew::ProcessRow","Data")
adb0e845	207	<<"Error in data given to the cluster finder"<<ENDLOG;
99a6d5c9	208	seq_mean = 1;
	209	seq_charge = 1;
	210	}
46fbeb8a	211
	212	//Calculate mean in pad direction:
	213	Int_t pad_mean = seq_charge*data[bin].fPad;
c3dd27a3	214	Int_t pad_error = data[bin].fPad*pad_mean;
b62fa469	215
46fbeb8a	216	//Compare with results on previous pad:
	217	for(UInt_t p=0; p<n_previous; p++)
	218	{
b62fa469	219	//dont merge sequences on the same pad twice
	220	if(previousPt[p]->fLastMergedPad==data[bin].fPad) continue;
	221
46fbeb8a	222	Int_t difference = seq_mean - previousPt[p]->fMean;
b62fa469	223	if(difference < -fMatch) break;
46fbeb8a	224
c3dd27a3	225	if(difference <= fMatch) //There is a match here!!
46fbeb8a	226	{
99a6d5c9	227	ClusterData *local = previousPt[p];
b62fa469	228
99a6d5c9	229	if(fDeconvPad)
	230	{
	231	if(seq_charge > local->fLastCharge)
	232	{
c3dd27a3	233	if(local->fChargeFalling) //The previous pad was falling
99a6d5c9	234	{
c3dd27a3	235	break; //create a new cluster
99a6d5c9	236	}
	237	}
	238	else
	239	local->fChargeFalling = 1;
	240	local->fLastCharge = seq_charge;
	241	}
	242
46fbeb8a	243	//Don't create a new cluster, because we found a match
	244	new_cluster = kFALSE;
	245
	246	//Update cluster on current pad with the matching one:
99a6d5c9	247	local->fTotalCharge += seq_charge;
99a6d5c9	248	local->fPad += pad_mean;
c3dd27a3	249	local->fPad2 += pad_error;
99a6d5c9	250	local->fTime += seq_average;
c3dd27a3	251	local->fTime2 += seq_error;
99a6d5c9	252	local->fMean = seq_mean;
99a6d5c9	253	local->fFlags++; //means we have more than one pad
b62fa469	254	local->fLastMergedPad = data[bin].fPad;
b62fa469	255
99a6d5c9	256	currentPt[n_current] = local;
46fbeb8a	257	n_current++;
46fbeb8a	258
46fbeb8a	259	break;
c3dd27a3	260	} //Checking for match at previous pad
c3dd27a3	261	} //Loop over results on previous pad.
46fbeb8a	262
	263	if(new_cluster)
	264	{
	265	//Start a new cluster. Add it to the clusterlist, and update
	266	//the list of pointers to clusters in current pad.
99a6d5c9	267	//current pad will be previous pad on next pad.
46fbeb8a	268
	269	//Add to the clusterlist:
	270	ClusterData *tmp = &clusterlist[n_total];
	271	tmp->fTotalCharge = seq_charge;
	272	tmp->fPad = pad_mean;
c3dd27a3	273	tmp->fPad2 = pad_error;
46fbeb8a	274	tmp->fTime = seq_average;
c3dd27a3	275	tmp->fTime2 = seq_error;
46fbeb8a	276	tmp->fMean = seq_mean;
46fbeb8a	277	tmp->fFlags = 0; //flags for 1 pad clusters
b62fa469	278	tmp->fLastMergedPad = data[bin].fPad;
99a6d5c9	279	if(fDeconvPad)
	280	{
	281	tmp->fChargeFalling = 0;
	282	tmp->fLastCharge = seq_charge;
	283	}
	284
46fbeb8a	285	//Update list of pointers to previous pad:
	286	currentPt[n_current] = &clusterlist[n_total];
	287	n_total++;
	288	n_current++;
	289	}
b62fa469	290
99a6d5c9	291	if(fDeconvTime)
99a6d5c9	292	if(new_bin >= 0) goto redo;
46fbeb8a	293	}//Loop over digits on this padrow
	294
	295	WriteClusters(n_total,clusterlist);
46fbeb8a	296	}
	297
	298	void AliL3ClustFinderNew::WriteClusters(Int_t n_clusters,ClusterData *list)
	299	{
	300	Int_t thisrow,thissector;
	301	UInt_t counter = fNClusters;
	302
	303	for(int j=0; j<n_clusters; j++)
	304	{
	305	if(!list[j].fFlags) continue; //discard 1 pad clusters
	306	if(list[j].fTotalCharge < fThreshold) continue; //noise cluster
d3b44290	307
d3b44290	308	Float_t xyz[3];
c3dd27a3	309	Float_t fpad =(Float_t)list[j].fPad /(Float_t)list[j].fTotalCharge;
	310	Float_t fpad2=fXYErr;
	311	Float_t ftime =(Float_t)list[j].fTime /(Float_t)list[j].fTotalCharge;
	312	Float_t ftime2=fZErr;
eeddc64d	313
c3dd27a3	314	if(fcalcerr) {
	315	fpad2=(Float_t)list[j].fPad2/(Float_t)list[j].fTotalCharge - fpad*fpad;
	316	fpad2 = sqrt(fpad2);
	317	ftime2=(Float_t)list[j].fTime2/(Float_t)list[j].fTotalCharge - ftime*ftime;
	318	ftime2 = sqrt(ftime2);
	319	}
	320
77f08cd3	321	if(fstdout==kTRUE)
b62fa469	322	cout<<"WriteCluster: padrow "<<fCurrentRow<<" pad "<<fpad << " +- "<<fpad2<<" time "<<ftime<<" +- "<<ftime2<<" charge "<<list[j].fTotalCharge<<endl;
d3b44290	323
494fad94	324	AliL3Transform::Slice2Sector(fCurrentSlice,fCurrentRow,thissector,thisrow);
494fad94	325	AliL3Transform::Raw2Local(xyz,thissector,thisrow,fpad,ftime);
46fbeb8a	326	if(xyz[0]==0) LOG(AliL3Log::kError,"AliL3ClustFinder","Cluster Finder")
	327	<<AliL3Log::kDec<<"Zero cluster"<<ENDLOG;
	328	if(fNClusters >= fMaxNClusters)
	329	{
	330	LOG(AliL3Log::kError,"AliL3ClustFinder::WriteClusters","Cluster Finder")
	331	<<AliL3Log::kDec<<"Too many clusters"<<ENDLOG;
	332	return;
d3b44290	333	}
95a00d93	334	fSpacePointData[counter].fCharge = list[j].fTotalCharge;
46fbeb8a	335	fSpacePointData[counter].fX = xyz[0];
	336	fSpacePointData[counter].fY = xyz[1];
	337	fSpacePointData[counter].fZ = xyz[2];
	338	fSpacePointData[counter].fPadRow = fCurrentRow;
c3dd27a3	339	fSpacePointData[counter].fXYErr = fpad2;
c3dd27a3	340	fSpacePointData[counter].fZErr = ftime2;
46fbeb8a	341	fSpacePointData[counter].fID = counter
d3b44290	342	+((fCurrentSlice&0x7f)<<25)+((fCurrentPatch&0x7)<<22);//Uli
10f815d9	343	#ifdef do_mc
	344	Int_t trackID[3];
	345	GetTrackID((Int_t)rint(fpad),(Int_t)rint(ftime),trackID);
6ab8a12a	346	//cout<<"padrow "<<fCurrentRow<<" pad "<<(Int_t)rint(fpad)<<" time "<<(Int_t)rint(ftime)<<" Trackid "<<trackID[0]<<endl;
10f815d9	347	fSpacePointData[counter].fTrackID[0] = trackID[0];
	348	fSpacePointData[counter].fTrackID[1] = trackID[1];
	349	fSpacePointData[counter].fTrackID[2] = trackID[2];
	350	#endif
9953e05f	351
46fbeb8a	352	fNClusters++;
46fbeb8a	353	counter++;
46fbeb8a	354	}
46fbeb8a	355	}
10f815d9	356
	357	#ifdef do_mc
	358	void AliL3ClustFinderNew::GetTrackID(Int_t pad,Int_t time,Int_t *trackID)
	359	{
	360	AliL3DigitRowData rowPt = (AliL3DigitRowData)fDigitRowData;
	361
355debe1	362	trackID[0]=trackID[1]=trackID[2]=-2;
10f815d9	363	//cout<<"Looking for pad "<<pad<<" time "<<time<<endl;
	364	for(Int_t i=fFirstRow; i<=fLastRow; i++)
	365	{
	366	if(rowPt->fRow < (UInt_t)fCurrentRow)
	367	{
	368	AliL3MemHandler::UpdateRowPointer(rowPt);
	369	continue;
	370	}
	371	AliL3DigitData digPt = (AliL3DigitData)rowPt->fDigitData;
	372	for(UInt_t j=0; j<rowPt->fNDigit; j++)
	373	{
	374	Int_t cpad = digPt[j].fPad;
	375	Int_t ctime = digPt[j].fTime;
355debe1	376	if(cpad != pad) continue;
	377	if(ctime != time) continue;
	378	//if(cpad != pad && ctime != ctime) continue;
10f815d9	379	//cout<<"Reading row "<<fCurrentRow<<" pad "<<cpad<<" time "<<ctime<<" trackID "<<digPt[j].fTrackID[0]<<endl;
	380	trackID[0] = digPt[j].fTrackID[0];
	381	trackID[1] = digPt[j].fTrackID[1];
	382	trackID[2] = digPt[j].fTrackID[2];
355debe1	383	break;
10f815d9	384	//cout<<"Reading trackID "<<trackID[0]<<endl;
	385	}
	386	break;
	387	}
	388
	389	}
	390	#endif