[u/mrichter/AliRoot.git] / HLT / misc / AliL3VHDLClusterFinder.cxx

// @(#) $Id$

// Author: Constantin Loizides <mailto:loizides@ikf.uni-frankfurt.de>
//*-- Copyright & copy ALICE HLT Group
/** \class AliL3VHDLClusterFinder
<pre>
//____________________________________________________
// AliL3VHDLClusterFinder
//
// The current VHDL cluster finder for HLT
// Based on STAR L3
//
// Most important parameters:
// fThreshold - threshold for noise clusters
// fMatch - length in time for overlapping sequences
</pre> 
*/

#include "AliL3StandardIncludes.h"
#include "AliL3RootTypes.h"
#include "AliL3Logging.h"
#include "AliL3AltroMemHandler.h"

//#define VHDLDEBUG
#include "AliL3VHDLClusterFinder.h"

#if __GNUC__ >= 3
using namespace std;
#endif


ClassImp(AliL3VHDLClusterFinder)

AliL3VHDLClusterFinder::AliL3VHDLClusterFinder()
{
  // default constructor
  fMatch = 4;
  fThreshold = 10;
  fMinMerge = 1;
  fNClusters=0;

  fXYErr = 0.2;
  fZErr = 0.3;
  fDeconvPad = kTRUE;
  fDeconvTime = kTRUE;
  fstdout = kFALSE;
  fcalcerr = kTRUE;

  Clear();
#ifdef VHDLDEBUG
  fdeb=fopen("vhdlclusterfinder.debug","w");
  //fdeb=stderr;
#endif
}

AliL3VHDLClusterFinder::~AliL3VHDLClusterFinder()
{
  // destructor
#ifdef VHDLDEBUG
  fclose(fdeb);
#endif
}

void AliL3VHDLClusterFinder::ProcessDigits()
{
  //Loop over data like the analyzer of the VHDL code
  const UChar_t kn=255;
  UShort_t rrow=0,rtime=0;
  UChar_t rpad=0,i=kn;
  UShort_t *charges=new UShort_t[kn];

  fNClusters=0;
  fRow=0;
  fPad=0;
  Clear();

  //loop over input data
  while(fAltromem.ReadSequence(rrow,rpad,rtime,i,&charges)){
#if 0
    cout << "Padrow " << (int)rrow << " pad " << (int)rpad << " time " <<(int)rtime << " charges ";
    for(UChar_t ii=0;ii<i;ii++) cout << (int)charges[ii] << " ";
    cout << endl;
#endif
#ifdef VHDLDEBUG
    fprintf(fdeb,"ProcessDigits: Input Data: %d %d %d charges:",(int)rrow,(int)rpad,(int)rtime);
    for(UChar_t ii=0;ii<i;ii++) fprintf(fdeb," %d",(int)charges[ii]);
    fprintf(fdeb,"\n");
#endif

    fNRow=rrow;
    fNPad=rpad;
    fTC=0,fMT=0,fST=0;

    //calculate sequence values 
    if(fDeconvTime){
      UChar_t ii=0;
      Int_t sl=0;
      Int_t charge=0,lcharge=0;
      Bool_t falling=kFALSE;
      while(ii<i){
	charge=charges[ii];
	if((falling)&&(charge>lcharge)){
	  fSM=rtime-sl/2;
	  MakeSequence();
	  ProcessSequence();
	  OutputMemory();

	  rtime=rtime-sl;      
	  falling=kFALSE;
	  sl=0;
	  fTC=0,fMT=0,fST=0;
	} else if(charge<lcharge) falling=kTRUE;

	fTC+=charge;
	fMT+=(rtime-sl)*charge;
	fST+=(rtime-sl)*(rtime-sl)*charge;
	sl++;
	ii++;
	lcharge=charge;
      } //end loop over sequence
      fSM=rtime-sl/2;
    } else { /* no deconvolution */
      for(UChar_t ii=0;ii<i;ii++){
	fTC+=charges[ii];
	fMT+=(rtime-ii)*charges[ii];
	fST+=(rtime-ii)*(rtime-ii)*charges[ii];
      }
      fSM=rtime-i/2;
    }

    MakeSequence();
    ProcessSequence();
    OutputMemory();      

#ifdef VHDLDEBUG
    fflush(fdeb);
#endif
    i=kn; //store size of charges array
  } //loop over data

  //flush everything left
  FlushMemory();
  while(fOP!=fFP) OutputMemory();
}

void AliL3VHDLClusterFinder::MakeSequence(){
  // makes the sequence
  if(!fTC) return;

  Int_t mp=fNPad*fTC;
  Int_t sp=fNPad*fNPad*fTC;

  fSeq.fTotalCharge=fTC;
  fSeq.fPad=mp;
  fSeq.fTime=fMT;
  fSeq.fPad2=sp;
  fSeq.fTime2=fST; 
  fSeq.fMean=0;
  fSeq.fMean=fSM;
  fSeq.fMerge=0;
  fSeq.fRow=fNRow;
  fSeq.fLastPad=fNPad;
  fSeq.fChargeFalling=0;
  if(fDeconvPad) fSeq.fLastCharge=fTC;
  else fSeq.fLastCharge=0;
}

void AliL3VHDLClusterFinder::ProcessSequence()
{
  // processes the sequence
  if(fNRow!=fRow) FlushMemory();
  else if(fNPad==fPad+1) PrepareMemory();
  else if(fNPad!=fPad) FlushMemory();

  //store new row and pad values
  fRow=fNRow;
  fPad=fNPad;

#ifdef VHDLDEBUG
  fprintf(fdeb,"ProcessSequence: Mean=%d TC=%d ",fSeq.fMean,fSeq.fTotalCharge);
#endif

  CompareSeq(); //merge or insert
}

void AliL3VHDLClusterFinder::PrepareMemory()
{
  // prepares the memory
#ifdef VHDLDEBUG
  fprintf(fdeb,"PrepareMemory %d %d %d\n",fRP,fEP,fWP);
#endif
  fRP=fEP;
  fEP=fWP;
}

void AliL3VHDLClusterFinder::FlushMemory()
{
  // flushes the memory
#ifdef VHDLDEBUG
  fprintf(fdeb,"FlushMemory %d %d %d %d\n",fFP,fRP,fEP,fWP);
#endif
  fFP=fWP;
  fRP=fWP;
  fEP=fWP;
}

void AliL3VHDLClusterFinder::CompareSeq()
{
  // compares sequences
  while(fRP!=fEP){
    Int_t diff=fSeqs[fPList[fRP]].fMean-fSeq.fMean;

    if(diff>fMatch){ //no match
      IncRPointer(); //cluster finished
      continue;
    } else if(diff<-fMatch){ //no match
      break;                 //insert new cluster       
    } else { //match found, merge it
      MergeSeq(); 
      return;
    }
  }

  InsertSeq(); //start new cluster
}

void AliL3VHDLClusterFinder::MergeSeq()
{
  // merges sequences
#ifdef VHDLDEBUG
  fprintf(fdeb,"merged with Mean=%d TC=%d (new Merge=%d) %d %d\n",fSeqs[fPList[fRP]].fMean,fSeqs[fPList[fRP]].fTotalCharge,fSeqs[fPList[fRP]].fMerge+1,fRow,fPad);
#endif
  if(fSeqs[fPList[fRP]].fRow!=fSeq.fRow){
    LOG(AliL3Log::kWarning,"AliL3VHDLClusterFinder::","Memory Check")
      <<"Sequences can be merged on the same rows only."<<ENDLOG;
  }
  if(fSeqs[fPList[fRP]].fLastPad+1!=fSeq.fLastPad){
    LOG(AliL3Log::kWarning,"AliL3VHDLClusterFinder::","Memory Check")
      <<"Sequences can be merged on consecutive pads only."<<ENDLOG;
  }
  
  if(fDeconvPad){
    if(fSeq.fTotalCharge > fSeqs[fPList[fRP]].fLastCharge){
      if(fSeqs[fPList[fRP]].fChargeFalling){ //The previous pad was falling
	IncRPointer();
	InsertSeq(); //start a new cluster
	return;
      }		    
    }
    else fSeqs[fPList[fRP]].fChargeFalling = 1;

    fSeqs[fPList[fRP]].fLastCharge = fSeq.fTotalCharge;
  }
  
  fSeqs[fPList[fRP]].fMean=fSeq.fMean; //take the new mean
  fSeqs[fPList[fRP]].fLastPad=fSeq.fLastPad;
  fSeqs[fPList[fRP]].fTotalCharge+=fSeq.fTotalCharge;
  fSeqs[fPList[fRP]].fPad+=fSeq.fPad;
  fSeqs[fPList[fRP]].fTime+=fSeq.fTime;
  fSeqs[fPList[fRP]].fPad2+=fSeq.fPad2;
  fSeqs[fPList[fRP]].fTime2+=fSeq.fTime2;
  fSeqs[fPList[fRP]].fMerge+=1;

  fPList[fWP]=fPList[fRP];
  fPList[fRP]=N_clmem; //mark it to be free

  IncRPointer();
  IncWPointer();
}

void AliL3VHDLClusterFinder::InsertSeq()
{
  // inserts sequence
#ifdef VHDLDEBUG
  fprintf(fdeb,"inserted %d %d\n",fRow,fPad);
#endif
  NextFreeIndex();    //get next index
  fSeqs[fFirst]=fSeq; //store data
  fPList[fWP]=fFirst; //store index
  IncWPointer();
}

void AliL3VHDLClusterFinder::OutputMemory()
{
  // output memory?
  Float_t mtime=0,mpad=0;
  Float_t mtime2=0,mpad2=0;
  UInt_t tc,row,mno;

  if(fOP!=fFP){
  //while(fOP!=fFP){

    UInt_t index=fPList[fOP];
    fPList[fOP]=N_clmem;
    //cout << fOP << " - " << fFP << " - " << index << endl;
    IncPointer(fOP);
    if(index>=N_clmem) return; //nothing to do
    //if(index>=N_clmem) continue; //nothing to do
    
    tc=fSeqs[index].fTotalCharge;
    mno=fSeqs[index].fMerge;
    row=fSeqs[index].fRow;
    if(tc!=0){
      mtime=(Float_t)(fSeqs[index].fTime)/tc;
      mtime2=sqrt((Float_t)(fSeqs[index].fTime2)/tc-mtime*mtime);
    }
    if(tc!=0){
      mpad=(Float_t)(fSeqs[index].fPad)/tc;
      mpad2=sqrt((Float_t)(fSeqs[index].fPad2)/tc-mpad*mpad);
    }

    if(mno<fMinMerge){ //noise cluster
#ifdef VHDLDEBUG
    fprintf(fdeb,"OutputMemory %d %d: noise cluster (merge): Row %d Pad %.3f Time %.3f TC %d Merge %d\n",fOP,fFP,row,mpad,mtime,tc,mno);
#endif          
      FreeSeq(index);  
      return;
      //continue;
    }

    if(tc<fThreshold){ //total charge below threshold
#ifdef VHDLDEBUG
    fprintf(fdeb,"OutputMemory %d %d: noise cluster (threshold): Row %d Pad %.3f Time %.3f TC %d Merge %d\n",fOP,fFP,row,mpad,mtime,tc,mno);
#endif          
      FreeSeq(index);
      return;
      //continue;
    }

#ifdef VHDLDEBUG
    fprintf(fdeb,"OutputMemory %d: Row %d Pad %.3f Time %.3f TC %d Merge %d\n",fNClusters,row,mpad,mtime,tc,mno);
#endif    

    //    if(stdout)
    //      cout<<"WriteCluster: padrow "<<row<<" pad "<<mpad<< " +- "<<mpad2<<" time "<<mtime<<" +- "<<mtime2<<" charge "<<tc<<endl;
    
    fNClusters++;
    FreeSeq(index);
  }
}

void AliL3VHDLClusterFinder::FreeSeq(UShort_t i)
{
  // frees the sequence
  ClearSeq(i);
  IncPointer(fLast,1,N_clmem);
}

void AliL3VHDLClusterFinder::Clear()
{
  // clears everything
  fFirst=0; //first list pointer
  fLast=0;  //last  list pointer

  fWP=0; //write pointer
  fRP=0; //read pointer
  fEP=0; //end read pointer
  fFP=0; //end flush pointer
  fOP=0; //output pointer

  fSeq.fRow=0;
  fSeq.fLastPad=0;
  fSeq.fTotalCharge=0;
  fSeq.fPad=0;
  fSeq.fTime=0;
  fSeq.fPad2=0;
  fSeq.fTime2=0; 
  fSeq.fMean=0;
  fSeq.fMerge=0;

  for(UInt_t i=0;i<N_mem;i++) fPList[i]=N_clmem;
  for(UInt_t i=0;i<N_clmem;i++) ClearSeq(i);
}

void AliL3VHDLClusterFinder::ClearSeq(UShort_t i){
  // clears a sequence
  fSeqs[i].fRow=0;
  fSeqs[i].fLastPad=0;
  fSeqs[i].fTotalCharge=0;
  fSeqs[i].fPad=0;
  fSeqs[i].fTime=0;
  fSeqs[i].fPad2=0;
  fSeqs[i].fTime2=0; 
  fSeqs[i].fMean=0;
  fSeqs[i].fMerge=0;
}

void AliL3VHDLClusterFinder::IncPointer(UShort_t &p, Short_t add, UShort_t N){
  // increments pointer by 'add' in a circular buffer 
  Short_t pp=p;
  pp+=add;  
  if(pp>=N) p=UShort_t(pp-N);
  else if(pp<0) p=UShort_t(pp+N);
  else p=UShort_t(pp);
}

void AliL3VHDLClusterFinder::IncRPointer(){
  // increments pointer fRP 
  IncPointer(fRP);
}

void AliL3VHDLClusterFinder::IncWPointer(){
  // increments pointer fWP
  IncPointer(fWP);

  if(fWP==fOP){
    LOG(AliL3Log::kWarning,"AliL3VHDLClusterFinder::IncWPointer","Memory Check")
      <<"Write pointer overwrites output pointer."<<ENDLOG;
  }
}

void AliL3VHDLClusterFinder::NextFreeIndex(){
  // finds next free index
  IncPointer(fFirst,1,N_clmem);
  if(fFirst==fLast) {
    LOG(AliL3Log::kFatal,"AliL3VHDLClusterFinder::GetFreeIndex","Memory Check")
      <<"No space left in sequence list: "<<fFirst<<"=="<<fLast<<ENDLOG;
  }
}

#if 0
void AliL3ClustFinderNew::WriteClusters(Int_t n_clusters,ClusterData *list)
{
  // writes clusters
  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

      fNClusters++;
      counter++;
    }
}
#endif
Commit	Line	Data
3e87ef69	1	// @(#) $Id$
1e843b89	2
3e87ef69	3	// Author: Constantin Loizides <mailto:loizides@ikf.uni-frankfurt.de>
3e87ef69	4	//*-- Copyright & copy ALICE HLT Group
8e1f5064	5	/** \class AliL3VHDLClusterFinder
3e87ef69	6	<pre>
8e1f5064	7	//____________________________________________________
	8	// AliL3VHDLClusterFinder
	9	//
	10	// The current VHDL cluster finder for HLT
	11	// Based on STAR L3
	12	//
	13	// Most important parameters:
	14	// fThreshold - threshold for noise clusters
	15	// fMatch - length in time for overlapping sequences
3e87ef69	16	</pre>
6aec85a4	17	*/
8e1f5064	18
54b54089	19	#include "AliL3StandardIncludes.h"
	20	#include "AliL3RootTypes.h"
	21	#include "AliL3Logging.h"
	22	#include "AliL3AltroMemHandler.h"
	23
	24	//#define VHDLDEBUG
	25	#include "AliL3VHDLClusterFinder.h"
	26
5929c18d	27	#if __GNUC__ >= 3
54b54089	28	using namespace std;
	29	#endif
	30
c8e2510a	31
8e1f5064	32	ClassImp(AliL3VHDLClusterFinder)
	33
	34	AliL3VHDLClusterFinder::AliL3VHDLClusterFinder()
	35	{
54b54089	36	// default constructor
8e1f5064	37	fMatch = 4;
	38	fThreshold = 10;
	39	fMinMerge = 1;
	40	fNClusters=0;
	41
	42	fXYErr = 0.2;
	43	fZErr = 0.3;
	44	fDeconvPad = kTRUE;
	45	fDeconvTime = kTRUE;
	46	fstdout = kFALSE;
	47	fcalcerr = kTRUE;
	48
	49	Clear();
3e87ef69	50	#ifdef VHDLDEBUG
8e1f5064	51	fdeb=fopen("vhdlclusterfinder.debug","w");
	52	//fdeb=stderr;
	53	#endif
	54	}
	55
	56	AliL3VHDLClusterFinder::~AliL3VHDLClusterFinder()
	57	{
54b54089	58	// destructor
3e87ef69	59	#ifdef VHDLDEBUG
8e1f5064	60	fclose(fdeb);
	61	#endif
	62	}
	63
	64	void AliL3VHDLClusterFinder::ProcessDigits()
	65	{
	66	//Loop over data like the analyzer of the VHDL code
54b54089	67	const UChar_t kn=255;
8e1f5064	68	UShort_t rrow=0,rtime=0;
54b54089	69	UChar_t rpad=0,i=kn;
54b54089	70	UShort_t *charges=new UShort_t[kn];
8e1f5064	71
	72	fNClusters=0;
	73	fRow=0;
	74	fPad=0;
	75	Clear();
	76
	77	//loop over input data
	78	while(fAltromem.ReadSequence(rrow,rpad,rtime,i,&charges)){
8e1f5064	79	#if 0
	80	cout << "Padrow " << (int)rrow << " pad " << (int)rpad << " time " <<(int)rtime << " charges ";
	81	for(UChar_t ii=0;ii<i;ii++) cout << (int)charges[ii] << " ";
	82	cout << endl;
	83	#endif
3e87ef69	84	#ifdef VHDLDEBUG
8e1f5064	85	fprintf(fdeb,"ProcessDigits: Input Data: %d %d %d charges:",(int)rrow,(int)rpad,(int)rtime);
	86	for(UChar_t ii=0;ii<i;ii++) fprintf(fdeb," %d",(int)charges[ii]);
	87	fprintf(fdeb,"\n");
	88	#endif
	89
	90	fNRow=rrow;
	91	fNPad=rpad;
c8e2510a	92	fTC=0,fMT=0,fST=0;
8e1f5064	93
8e1f5064	94	//calculate sequence values
c8e2510a	95	if(fDeconvTime){
	96	UChar_t ii=0;
	97	Int_t sl=0;
	98	Int_t charge=0,lcharge=0;
	99	Bool_t falling=kFALSE;
	100	while(ii<i){
	101	charge=charges[ii];
	102	if((falling)&&(charge>lcharge)){
	103	fSM=rtime-sl/2;
	104	MakeSequence();
	105	ProcessSequence();
	106	OutputMemory();
	107
	108	rtime=rtime-sl;
	109	falling=kFALSE;
	110	sl=0;
	111	fTC=0,fMT=0,fST=0;
	112	} else if(charge<lcharge) falling=kTRUE;
	113
	114	fTC+=charge;
	115	fMT+=(rtime-sl)*charge;
	116	fST+=(rtime-sl)(rtime-sl)charge;
	117	sl++;
	118	ii++;
	119	lcharge=charge;
	120	} //end loop over sequence
	121	fSM=rtime-sl/2;
	122	} else { /* no deconvolution */
	123	for(UChar_t ii=0;ii<i;ii++){
	124	fTC+=charges[ii];
	125	fMT+=(rtime-ii)*charges[ii];
	126	fST+=(rtime-ii)(rtime-ii)charges[ii];
	127	}
	128	fSM=rtime-i/2;
8e1f5064	129	}
c8e2510a	130
c8e2510a	131	MakeSequence();
8e1f5064	132	ProcessSequence();
c8e2510a	133	OutputMemory();
8e1f5064	134
3e87ef69	135	#ifdef VHDLDEBUG
8e1f5064	136	fflush(fdeb);
8e1f5064	137	#endif
54b54089	138	i=kn; //store size of charges array
8e1f5064	139	} //loop over data
	140
	141	//flush everything left
	142	FlushMemory();
	143	while(fOP!=fFP) OutputMemory();
	144	}
	145
c8e2510a	146	void AliL3VHDLClusterFinder::MakeSequence(){
54b54089	147	// makes the sequence
c8e2510a	148	if(!fTC) return;
	149
	150	Int_t mp=fNPad*fTC;
	151	Int_t sp=fNPadfNPadfTC;
	152
	153	fSeq.fTotalCharge=fTC;
	154	fSeq.fPad=mp;
	155	fSeq.fTime=fMT;
	156	fSeq.fPad2=sp;
	157	fSeq.fTime2=fST;
	158	fSeq.fMean=0;
	159	fSeq.fMean=fSM;
	160	fSeq.fMerge=0;
	161	fSeq.fRow=fNRow;
	162	fSeq.fLastPad=fNPad;
	163	fSeq.fChargeFalling=0;
	164	if(fDeconvPad) fSeq.fLastCharge=fTC;
	165	else fSeq.fLastCharge=0;
	166	}
	167
8e1f5064	168	void AliL3VHDLClusterFinder::ProcessSequence()
8e1f5064	169	{
54b54089	170	// processes the sequence
8e1f5064	171	if(fNRow!=fRow) FlushMemory();
	172	else if(fNPad==fPad+1) PrepareMemory();
	173	else if(fNPad!=fPad) FlushMemory();
	174
	175	//store new row and pad values
	176	fRow=fNRow;
	177	fPad=fNPad;
	178
3e87ef69	179	#ifdef VHDLDEBUG
8e1f5064	180	fprintf(fdeb,"ProcessSequence: Mean=%d TC=%d ",fSeq.fMean,fSeq.fTotalCharge);
	181	#endif
	182
	183	CompareSeq(); //merge or insert
	184	}
	185
	186	void AliL3VHDLClusterFinder::PrepareMemory()
	187	{
54b54089	188	// prepares the memory
3e87ef69	189	#ifdef VHDLDEBUG
8e1f5064	190	fprintf(fdeb,"PrepareMemory %d %d %d\n",fRP,fEP,fWP);
	191	#endif
	192	fRP=fEP;
	193	fEP=fWP;
	194	}
	195
	196	void AliL3VHDLClusterFinder::FlushMemory()
	197	{
54b54089	198	// flushes the memory
3e87ef69	199	#ifdef VHDLDEBUG
8e1f5064	200	fprintf(fdeb,"FlushMemory %d %d %d %d\n",fFP,fRP,fEP,fWP);
	201	#endif
	202	fFP=fWP;
	203	fRP=fWP;
	204	fEP=fWP;
	205	}
	206
	207	void AliL3VHDLClusterFinder::CompareSeq()
	208	{
54b54089	209	// compares sequences
8e1f5064	210	while(fRP!=fEP){
	211	Int_t diff=fSeqs[fPList[fRP]].fMean-fSeq.fMean;
	212
	213	if(diff>fMatch){ //no match
	214	IncRPointer(); //cluster finished
	215	continue;
	216	} else if(diff<-fMatch){ //no match
	217	break; //insert new cluster
c8e2510a	218	} else { //match found, merge it
8e1f5064	219	MergeSeq();
	220	return;
	221	}
	222	}
	223
	224	InsertSeq(); //start new cluster
	225	}
	226
	227	void AliL3VHDLClusterFinder::MergeSeq()
	228	{
54b54089	229	// merges sequences
3e87ef69	230	#ifdef VHDLDEBUG
6aec85a4	231	fprintf(fdeb,"merged with Mean=%d TC=%d (new Merge=%d) %d %d\n",fSeqs[fPList[fRP]].fMean,fSeqs[fPList[fRP]].fTotalCharge,fSeqs[fPList[fRP]].fMerge+1,fRow,fPad);
8e1f5064	232	#endif
fbd800c2	233	if(fSeqs[fPList[fRP]].fRow!=fSeq.fRow){
8e1f5064	234	LOG(AliL3Log::kWarning,"AliL3VHDLClusterFinder::","Memory Check")
	235	<<"Sequences can be merged on the same rows only."<<ENDLOG;
	236	}
	237	if(fSeqs[fPList[fRP]].fLastPad+1!=fSeq.fLastPad){
	238	LOG(AliL3Log::kWarning,"AliL3VHDLClusterFinder::","Memory Check")
	239	<<"Sequences can be merged on consecutive pads only."<<ENDLOG;
	240	}
c8e2510a	241
	242	if(fDeconvPad){
	243	if(fSeq.fTotalCharge > fSeqs[fPList[fRP]].fLastCharge){
	244	if(fSeqs[fPList[fRP]].fChargeFalling){ //The previous pad was falling
	245	IncRPointer();
	246	InsertSeq(); //start a new cluster
	247	return;
	248	}
	249	}
	250	else fSeqs[fPList[fRP]].fChargeFalling = 1;
8e1f5064	251
c8e2510a	252	fSeqs[fPList[fRP]].fLastCharge = fSeq.fTotalCharge;
	253	}
	254
8e1f5064	255	fSeqs[fPList[fRP]].fMean=fSeq.fMean; //take the new mean
	256	fSeqs[fPList[fRP]].fLastPad=fSeq.fLastPad;
	257	fSeqs[fPList[fRP]].fTotalCharge+=fSeq.fTotalCharge;
	258	fSeqs[fPList[fRP]].fPad+=fSeq.fPad;
	259	fSeqs[fPList[fRP]].fTime+=fSeq.fTime;
	260	fSeqs[fPList[fRP]].fPad2+=fSeq.fPad2;
	261	fSeqs[fPList[fRP]].fTime2+=fSeq.fTime2;
	262	fSeqs[fPList[fRP]].fMerge+=1;
	263
	264	fPList[fWP]=fPList[fRP];
	265	fPList[fRP]=N_clmem; //mark it to be free
	266
	267	IncRPointer();
	268	IncWPointer();
	269	}
	270
	271	void AliL3VHDLClusterFinder::InsertSeq()
	272	{
54b54089	273	// inserts sequence
3e87ef69	274	#ifdef VHDLDEBUG
6aec85a4	275	fprintf(fdeb,"inserted %d %d\n",fRow,fPad);
8e1f5064	276	#endif
	277	NextFreeIndex(); //get next index
	278	fSeqs[fFirst]=fSeq; //store data
	279	fPList[fWP]=fFirst; //store index
	280	IncWPointer();
	281	}
	282
	283	void AliL3VHDLClusterFinder::OutputMemory()
	284	{
54b54089	285	// output memory?
8e1f5064	286	Float_t mtime=0,mpad=0;
	287	Float_t mtime2=0,mpad2=0;
	288	UInt_t tc,row,mno;
	289
	290	if(fOP!=fFP){
	291	//while(fOP!=fFP){
	292
	293	UInt_t index=fPList[fOP];
	294	fPList[fOP]=N_clmem;
	295	//cout << fOP << " - " << fFP << " - " << index << endl;
	296	IncPointer(fOP);
	297	if(index>=N_clmem) return; //nothing to do
	298	//if(index>=N_clmem) continue; //nothing to do
	299
	300	tc=fSeqs[index].fTotalCharge;
	301	mno=fSeqs[index].fMerge;
	302	row=fSeqs[index].fRow;
	303	if(tc!=0){
	304	mtime=(Float_t)(fSeqs[index].fTime)/tc;
	305	mtime2=sqrt((Float_t)(fSeqs[index].fTime2)/tc-mtime*mtime);
	306	}
	307	if(tc!=0){
	308	mpad=(Float_t)(fSeqs[index].fPad)/tc;
	309	mpad2=sqrt((Float_t)(fSeqs[index].fPad2)/tc-mpad*mpad);
	310	}
	311
	312	if(mno<fMinMerge){ //noise cluster
3e87ef69	313	#ifdef VHDLDEBUG
8e1f5064	314	fprintf(fdeb,"OutputMemory %d %d: noise cluster (merge): Row %d Pad %.3f Time %.3f TC %d Merge %d\n",fOP,fFP,row,mpad,mtime,tc,mno);
	315	#endif
	316	FreeSeq(index);
	317	return;
	318	//continue;
	319	}
	320
	321	if(tc<fThreshold){ //total charge below threshold
3e87ef69	322	#ifdef VHDLDEBUG
8e1f5064	323	fprintf(fdeb,"OutputMemory %d %d: noise cluster (threshold): Row %d Pad %.3f Time %.3f TC %d Merge %d\n",fOP,fFP,row,mpad,mtime,tc,mno);
	324	#endif
	325	FreeSeq(index);
	326	return;
	327	//continue;
	328	}
	329
3e87ef69	330	#ifdef VHDLDEBUG
8e1f5064	331	fprintf(fdeb,"OutputMemory %d: Row %d Pad %.3f Time %.3f TC %d Merge %d\n",fNClusters,row,mpad,mtime,tc,mno);
	332	#endif
	333
e33f3609	334	// if(stdout)
e33f3609	335	// cout<<"WriteCluster: padrow "<<row<<" pad "<<mpad<< " +- "<<mpad2<<" time "<<mtime<<" +- "<<mtime2<<" charge "<<tc<<endl;
8e1f5064	336
	337	fNClusters++;
	338	FreeSeq(index);
	339	}
	340	}
	341
62bb4b3d	342	void AliL3VHDLClusterFinder::FreeSeq(UShort_t i)
8e1f5064	343	{
54b54089	344	// frees the sequence
8e1f5064	345	ClearSeq(i);
	346	IncPointer(fLast,1,N_clmem);
	347	}
	348
	349	void AliL3VHDLClusterFinder::Clear()
	350	{
54b54089	351	// clears everything
8e1f5064	352	fFirst=0; //first list pointer
	353	fLast=0; //last list pointer
	354
	355	fWP=0; //write pointer
	356	fRP=0; //read pointer
	357	fEP=0; //end read pointer
	358	fFP=0; //end flush pointer
	359	fOP=0; //output pointer
	360
	361	fSeq.fRow=0;
	362	fSeq.fLastPad=0;
	363	fSeq.fTotalCharge=0;
	364	fSeq.fPad=0;
	365	fSeq.fTime=0;
	366	fSeq.fPad2=0;
	367	fSeq.fTime2=0;
	368	fSeq.fMean=0;
	369	fSeq.fMerge=0;
	370
	371	for(UInt_t i=0;i<N_mem;i++) fPList[i]=N_clmem;
	372	for(UInt_t i=0;i<N_clmem;i++) ClearSeq(i);
	373	}
	374
	375	void AliL3VHDLClusterFinder::ClearSeq(UShort_t i){
54b54089	376	// clears a sequence
8e1f5064	377	fSeqs[i].fRow=0;
	378	fSeqs[i].fLastPad=0;
	379	fSeqs[i].fTotalCharge=0;
	380	fSeqs[i].fPad=0;
	381	fSeqs[i].fTime=0;
	382	fSeqs[i].fPad2=0;
	383	fSeqs[i].fTime2=0;
	384	fSeqs[i].fMean=0;
	385	fSeqs[i].fMerge=0;
	386	}
	387
	388	void AliL3VHDLClusterFinder::IncPointer(UShort_t &p, Short_t add, UShort_t N){
54b54089	389	// increments pointer by 'add' in a circular buffer
8e1f5064	390	Short_t pp=p;
	391	pp+=add;
	392	if(pp>=N) p=UShort_t(pp-N);
	393	else if(pp<0) p=UShort_t(pp+N);
	394	else p=UShort_t(pp);
	395	}
	396
62bb4b3d	397	void AliL3VHDLClusterFinder::IncRPointer(){
54b54089	398	// increments pointer fRP
8e1f5064	399	IncPointer(fRP);
	400	}
	401
62bb4b3d	402	void AliL3VHDLClusterFinder::IncWPointer(){
54b54089	403	// increments pointer fWP
8e1f5064	404	IncPointer(fWP);
	405
	406	if(fWP==fOP){
	407	LOG(AliL3Log::kWarning,"AliL3VHDLClusterFinder::IncWPointer","Memory Check")
	408	<<"Write pointer overwrites output pointer."<<ENDLOG;
	409	}
	410	}
	411
62bb4b3d	412	void AliL3VHDLClusterFinder::NextFreeIndex(){
54b54089	413	// finds next free index
8e1f5064	414	IncPointer(fFirst,1,N_clmem);
	415	if(fFirst==fLast) {
	416	LOG(AliL3Log::kFatal,"AliL3VHDLClusterFinder::GetFreeIndex","Memory Check")
	417	<<"No space left in sequence list: "<<fFirst<<"=="<<fLast<<ENDLOG;
	418	}
	419	}
	420
	421	#if 0
	422	void AliL3ClustFinderNew::WriteClusters(Int_t n_clusters,ClusterData *list)
	423	{
54b54089	424	// writes clusters
8e1f5064	425	Int_t thisrow,thissector;
	426	UInt_t counter = fNClusters;
	427
	428	for(int j=0; j<n_clusters; j++)
	429	{
	430	if(!list[j].fFlags) continue; //discard 1 pad clusters
	431	if(list[j].fTotalCharge < fThreshold) continue; //noise cluster
	432
	433	Float_t xyz[3];
	434	Float_t fpad =(Float_t)list[j].fPad /(Float_t)list[j].fTotalCharge;
	435	Float_t fpad2=fXYErr;
	436	Float_t ftime =(Float_t)list[j].fTime /(Float_t)list[j].fTotalCharge;
	437	Float_t ftime2=fZErr;
	438
	439	if(fcalcerr) {
	440	fpad2=(Float_t)list[j].fPad2/(Float_t)list[j].fTotalCharge - fpad*fpad;
	441	fpad2 = sqrt(fpad2);
	442	ftime2=(Float_t)list[j].fTime2/(Float_t)list[j].fTotalCharge - ftime*ftime;
	443	ftime2 = sqrt(ftime2);
	444	}
	445
	446	if(fstdout==kTRUE)
	447	cout<<"WriteCluster: padrow "<<fCurrentRow<<" pad "<<fpad << "+-"<<fpad2<<" time "<<ftime<<"+-"<<ftime2<<" charge "<<list[j].fTotalCharge<<endl;
	448
	449	AliL3Transform::Slice2Sector(fCurrentSlice,fCurrentRow,thissector,thisrow);
	450	AliL3Transform::Raw2Local(xyz,thissector,thisrow,fpad,ftime);
	451	if(xyz[0]==0) LOG(AliL3Log::kError,"AliL3ClustFinder","Cluster Finder")
	452	<<AliL3Log::kDec<<"Zero cluster"<<ENDLOG;
	453	if(fNClusters >= fMaxNClusters)
	454	{
	455	LOG(AliL3Log::kError,"AliL3ClustFinder::WriteClusters","Cluster Finder")
	456	<<AliL3Log::kDec<<"Too many clusters"<<ENDLOG;
	457	return;
	458	}
	459	fSpacePointData[counter].fCharge = list[j].fTotalCharge;
	460	fSpacePointData[counter].fX = xyz[0];
	461	fSpacePointData[counter].fY = xyz[1];
	462	fSpacePointData[counter].fZ = xyz[2];
	463	fSpacePointData[counter].fPadRow = fCurrentRow;
	464	fSpacePointData[counter].fXYErr = fpad2;
	465	fSpacePointData[counter].fZErr = ftime2;
	466	fSpacePointData[counter].fID = counter
	467	+((fCurrentSlice&0x7f)<<25)+((fCurrentPatch&0x7)<<22);//Uli
	468
	469	fNClusters++;
	470	counter++;
	471	}
	472	}
	473	#endif
	474