/************************************************************************** * This file is property of and copyright by the ALICE HLT Project * * All rights reserved. * * * * Primary Authors: * * Indranil Das * * * * 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$ */ /////////////////////////////////////////////// //Author : Indranil Das, SINP, INDIA // Sukalyan Chattopadhyay, SINP, INDIA // //Email : indra.das@saha.ac.in // sukalyan.chattopadhyay@saha.ac.in // // This class implements a hit reconstruction algorithm for the dimuon // high level trigger. // The algorithm finds 3 pad clusters by looking for unique pads with a charge // above a certain threshold. A centre of gravity type calculation is applied // to the three pads forming the cluster to find the hit's X or Y coordinate // along the non-bending and bending planes individually. // The sepperate X and Y coordinates are then merged to give the full coordinate // of the hit. ///////////////////////////////////////////////// #include "AliHLTMUONHitReconstructor.h" #include "AliHLTMUONRecHitsBlockStruct.h" #include const int AliHLTMUONHitReconstructor::fgkDetectorId = 0xA00; const int AliHLTMUONHitReconstructor::fgkDDLOffSet = 12 ; const int AliHLTMUONHitReconstructor::fgkNofDDL = 8 ; const int AliHLTMUONHitReconstructor::fgkDDLHeaderSize = 8; const int AliHLTMUONHitReconstructor::fgkEvenLutSize = 1645632 + 1; const int AliHLTMUONHitReconstructor::fgkOddLutSize = 3363840 + 1; const int AliHLTMUONHitReconstructor::fgkLutLine[2] = {54208, 59648}; const int AliHLTMUONHitReconstructor::fgkMinIdManuChannel[2] = {917696, 64}; const int AliHLTMUONHitReconstructor::fgkMaxIdManuChannel[2] = {2563327, 3363903}; const float AliHLTMUONHitReconstructor::fgkHalfPadSize[3] = {1.25, 2.50, 5.00}; AliHLTMUONHitReconstructor::AliHLTMUONHitReconstructor(): fkBlockHeaderSize(8), fkDspHeaderSize(8), fkBuspatchHeaderSize(4), fDCCut(0), fPadData(NULL), fLookUpTableData(NULL), fRecPoints(NULL), fRecPointsCount(NULL), fMaxRecPointsCount(0), fCentralCountB(0), fCentralCountNB(0), fIdOffSet(0), fDDLId(0), fDigitPerDDL(0), fDetManuChannelIdList(NULL), fCentralChargeB(NULL), fCentralChargeNB(NULL), fRecX(NULL), fRecY(NULL), fAvgChargeX(NULL), fAvgChargeY(NULL), fNofFiredDetElem(0), fDebugLevel(0), fBusToDetElem(), fBusToDDL() { // ctor if(AliHLTMUONHitReconstructor::fgkEvenLutSize > AliHLTMUONHitReconstructor::fgkOddLutSize){ fPadData = new DHLTPad[AliHLTMUONHitReconstructor::fgkEvenLutSize]; } else{ fPadData = new DHLTPad[AliHLTMUONHitReconstructor::fgkOddLutSize]; } fkBlockHeaderSize = 8; fkDspHeaderSize = 8; fkBuspatchHeaderSize = 4; bzero(fGetIdTotalData,336*80*2*sizeof(int)); } AliHLTMUONHitReconstructor::~AliHLTMUONHitReconstructor() { // dtor if(fPadData){ delete []fPadData; fPadData = NULL; } if(fLookUpTableData){ delete []fLookUpTableData; fLookUpTableData = NULL; } } int AliHLTMUONHitReconstructor::GetLutLine(int iDDL) const { return ( iDDL<16 ) ? fgkLutLine[0] : fgkLutLine[1] ; } bool AliHLTMUONHitReconstructor::LoadLookUpTable(DHLTLut* lookUpTableData, int lookUpTableId) { // function that loads LookUpTable (= position of each pad with electronic channel associated with it) if(lookUpTableId= fgkDDLOffSet + fgkNofDDL){ HLTError("DDL number is out of range (must be %d<=iDDL<%d)\n",fgkDDLOffSet,fgkDDLOffSet+fgkNofDDL); return false; } fDDLId = lookUpTableId; int lutSize = ((lookUpTableId%2)==0) ? fgkEvenLutSize : fgkOddLutSize ; int nofLutLine = GetLutLine(lookUpTableId); int idOffSet = fgkMinIdManuChannel[lookUpTableId%2]; int detManuChannelId; fLookUpTableData = new DHLTLut[lutSize]; fLookUpTableData[0].fIdManuChannel = 0; fLookUpTableData[0].fIX = 0 ; fLookUpTableData[0].fIY = 0 ; fLookUpTableData[0].fRealX = 0.0 ; fLookUpTableData[0].fRealY = 0.0 ; fLookUpTableData[0].fRealZ = 0.0 ; fLookUpTableData[0].fPlane = -1 ; fLookUpTableData[0].fPcbZone = -1 ; for(int i=0; i(rawData); const AliHLTUInt32_t* buffer = rawData; fIdOffSet= fgkMinIdManuChannel[(fDDLId%2)]; fDetManuChannelIdList = new int[rawDataSize]; int index = 0; int dataCount = 0; fNofFiredDetElem = 0; int detElemId = 0 ; int buspatchId = 0; int prevDetElemId = 0 ; int totalBlockSize,blockRawDataSize; int totalDspSize,dspRawDataSize; int totalBuspatchSize,buspatchRawDataSize; int indexDsp,indexBuspatch,indexRawData; unsigned int dataWord; int charge; int idManuChannel; for(int iBlock = 0; iBlock < 2 ;iBlock++){ // loop over 2 blocks totalBlockSize = buffer[index + 1]; blockRawDataSize = buffer[index + 2]; indexDsp = index + fkBlockHeaderSize; HLTDebug("totalBlockSize : %d, blockRawDataSize : %d",totalBlockSize,blockRawDataSize); while(blockRawDataSize > 0){ totalDspSize = buffer[indexDsp + 1]; dspRawDataSize = buffer[indexDsp + 2]; HLTDebug(" totalDSPSize : %d, dspRawDataSize : %d",totalDspSize,dspRawDataSize); //if(buffer[indexDsp+1] == 1) dspRawDataSize --; // temporary solution to read buspatches indexBuspatch = indexDsp + fkDspHeaderSize + 2; // this extra 2 word comes from the faulty defination of Dsp header size while(dspRawDataSize > 0){ totalBuspatchSize = buffer[indexBuspatch + 1]; buspatchRawDataSize = buffer[indexBuspatch + 2]; buspatchId = buffer[indexBuspatch + 3]; if((detElemId = fBusToDetElem[buspatchId])==0){ HLTError("No Detection element found for buspatch : %d",buspatchId); return false; } HLTDebug("\ttotalBusPatchSize : %d, buspatchRawDataSize : %d",totalBuspatchSize,buspatchRawDataSize); indexRawData = indexBuspatch + fkBuspatchHeaderSize; while(buspatchRawDataSize > 0){ dataWord = buffer[indexRawData]; charge = (unsigned short)(dataWord & 0xFFF); idManuChannel = 0x0; idManuChannel = (idManuChannel|(detElemId%100))<<17; idManuChannel |= (dataWord >> 12) & 0x1FFFF; idManuChannel -= fIdOffSet ; if(charge > fDCCut && charge > 5){ // (charge > 4) is due cut out the noise level fPadData[idManuChannel].fBuspatchId = buspatchId; fPadData[idManuChannel].fDetElemId = detElemId; fPadData[idManuChannel].fIdManuChannel = idManuChannel; fPadData[idManuChannel].fIX = fLookUpTableData[idManuChannel+1].fIX; fPadData[idManuChannel].fIY = fLookUpTableData[idManuChannel+1].fIY; fPadData[idManuChannel].fRealX = fLookUpTableData[idManuChannel+1].fRealX; fPadData[idManuChannel].fRealY = fLookUpTableData[idManuChannel+1].fRealY; fPadData[idManuChannel].fRealZ = fLookUpTableData[idManuChannel+1].fRealZ; fPadData[idManuChannel].fPcbZone = fLookUpTableData[idManuChannel+1].fPcbZone; fPadData[idManuChannel].fPlane = fLookUpTableData[idManuChannel+1].fPlane; fPadData[idManuChannel].fCharge = charge; fDetManuChannelIdList[dataCount] = idManuChannel; if(detElemId != prevDetElemId){ if(fNofFiredDetElem>0){ fMaxFiredPerDetElem[fNofFiredDetElem-1] = dataCount; } fNofFiredDetElem++; prevDetElemId = detElemId ; } HLTDebug("\t buspatch : %d, detele : %d, id : %d, manu : %d, channel : %d, X : %f, Y: %f", fPadData[idManuChannel].fBuspatchId,fPadData[idManuChannel].fDetElemId, idManuChannel,((dataWord >> 12) & 0x7FF),((dataWord >> 23) & 0x3F), fPadData[idManuChannel].fRealX,fPadData[idManuChannel].fRealY); dataCount ++; }// if charge is more than DC Cut limit condition indexRawData++; buspatchRawDataSize --; } indexBuspatch += totalBuspatchSize; dspRawDataSize -= totalBuspatchSize; }// buspatch loop indexDsp += totalDspSize; blockRawDataSize -= totalDspSize; }// DSP loop index = totalBlockSize; }// Block loop fDigitPerDDL = dataCount; fMaxFiredPerDetElem[fNofFiredDetElem-1] = dataCount; if(fDigitPerDDL == 0){ HLTInfo("An Empty DDL File found"); } return true; } bool AliHLTMUONHitReconstructor::FindRecHits() { // fuction that calls hit reconstruction detector element-wise for(int iDet=0; iDet0) FindCentralHits(fMaxFiredPerDetElem[iDet-1],fMaxFiredPerDetElem[iDet]); else FindCentralHits(0,fMaxFiredPerDetElem[iDet]); RecXRecY(); if(!MergeRecHits()){ HLTError("Failed to merge hits\n"); return false; } if(iDet==0) for(int i=0;i (0+1) && fPadData[idManuChannel].fIY < (79 - 1)){ //if(fPadData[idManuChannel].fIY > 0){ if(fCentralCountB>0){ hasFind = false; for(b = 0;b fPadData[idManuChannelCentral].fCharge){ fCentralChargeB[b] = idManuChannel; }// if condn on pad charge }// if condon on pad position }// for loop over b if(!hasFind){ fCentralChargeB[fCentralCountB] = idManuChannel; fCentralCountB++; } } else{ fCentralChargeB[fCentralCountB] = idManuChannel; fCentralCountB++; }// check the size of centralHitB for(b = 0;b 2 (to avoid edge value pb) }// B/Nb checking else{ if(fCentralCountNB>0){ hasFind = false; for(nb = 0;nb fPadData[idManuChannelCentral].fCharge){ fCentralChargeNB[nb] = idManuChannel; }// if condn over to find higher charge }// if condn over to find position }// for loop over presently all nb values if(!hasFind){ fCentralChargeNB[fCentralCountNB] = idManuChannel; fCentralCountNB++; } }// centralHitNB size test else{ fCentralChargeNB[fCentralCountNB] = idManuChannel; fCentralCountNB++; }// centralHitNB size test }// fill for bending and nonbending hit }// detElemId loop } void AliHLTMUONHitReconstructor::RecXRecY() { // find reconstructed X and Y for each plane separately int b,nb; int idCentral; int idLower = 0; int idUpper = 0; int idRight = 0; int idLeft = 0; fRecY = new float[fCentralCountB]; fRecX = new float[fCentralCountNB]; fAvgChargeY = new float[fCentralCountB]; fAvgChargeX = new float[fCentralCountNB]; for(b=0;b fAvgChargeY[j]){ fRecY[i] = (fRecY[i]*fAvgChargeY[i] + fRecY[j]*fAvgChargeY[j] )/(fAvgChargeY[i] + fAvgChargeY[j]); fRecY[j] = 0.0; } else{ fRecY[j] = (fRecY[i]*fAvgChargeY[i] + fRecY[j]*fAvgChargeY[j] )/(fAvgChargeY[i] + fAvgChargeY[j]); fRecY[i] = 0.0; }// search for higher charge }//pad position }//j for loop }//if fRecY[i] != 0.0 }// i for loop // MERGE Non Bending Plane hits, which are place side by side for(int i=0;i fAvgChargeX[j]){ fRecX[i] = (fRecX[i]*fAvgChargeX[i] + fRecX[j]*fAvgChargeX[j] )/(fAvgChargeX[i] + fAvgChargeX[j]); fRecX[j] = 0.0; } else{ fRecX[j] = (fRecX[i]*fAvgChargeX[i] + fRecX[j]*fAvgChargeX[j] )/(fAvgChargeX[i] + fAvgChargeX[j]); fRecX[i] = 0.0; }// search for higher charge }//pad position }//j for loop }//if fRecX[i] != 0.0 }// i for loop // Merge bending Plane hits with Non Bending for(int b=0;b fabsf(padCenterXB)) diffX = fabsf(fRecX[nb]) - fabsf(padCenterXB); else diffX = fabsf(padCenterXB) - fabsf(fRecX[nb]); if(fabsf(padCenterYNB)>fabsf(fRecY[b])) diffY = fabsf(padCenterYNB) - fabsf(fRecY[b]); else diffY = fabsf(fRecY[b]) - fabsf(padCenterYNB); if(diffX < halfPadLengthX && diffY < halfPadLengthY ){//&& fPadData[idCentralB].fIY != 0){ // First check that we have not overflowed the buffer. if((*fRecPointsCount) == fMaxRecPointsCount){ HLTError("Nof RecHit (i.e. %d) exceeds the max nof RecHit limit %d\n",(*fRecPointsCount),fMaxRecPointsCount); return false; } //fRecPoints[(*fRecPointsCount)].fId = idCentralB; fRecPoints[(*fRecPointsCount)].fX = fRecX[nb]; fRecPoints[(*fRecPointsCount)].fY = fRecY[b]; fRecPoints[(*fRecPointsCount)].fZ = fPadData[idCentralB].fRealZ; //fRecPoints[(*fRecPointsCount)].fDetElemId = (AliHLTUInt32_t)fPadData[idCentralB].fDetElemId; (*fRecPointsCount)++; }//if lies wihtin 5.0 mm }// condn over fRecX ! = 0.0 }// loop over NB side }// condn on fRecY[b] != 0.0 }// loop over B side; if(fRecX){ delete []fRecX; fRecX = NULL; } if(fRecY){ delete []fRecY; fRecY = NULL; } if(fAvgChargeX){ delete []fAvgChargeX; fAvgChargeX = NULL; } if(fAvgChargeY){ delete []fAvgChargeY; fAvgChargeY = NULL; } return true; } void AliHLTMUONHitReconstructor::Clear() { for(int iPad=0;iPad