HLT/MUON/OnlineAnalysis/AliHLTMUONHitReconstructor.cxx

   1 /**************************************************************************
   2  * This file is property of and copyright by the ALICE HLT Project        *
   3  * All rights reserved.                                                   *
   4  *                                                                        *
   5  * Primary Authors:                                                       *
   6  *   Indranil Das <indra.das@saha.ac.in>                                  *
   7  *                                                                        *
   8  * Permission to use, copy, modify and distribute this software and its   *
   9  * documentation strictly for non-commercial purposes is hereby granted   *
  10  * without fee, provided that the above copyright notice appears in all   *
  11  * copies and that both the copyright notice and this permission notice   *
  12  * appear in the supporting documentation. The authors make no claims     *
  13  * about the suitability of this software for any purpose. It is          *
  14  * provided "as is" without express or implied warranty.                  *
  15  **************************************************************************/
  16
  17 /* $Id$ */
  18
  19 ///////////////////////////////////////////////
  20 //Author : Indranil Das, SINP, INDIA
  21 //         Sukalyan Chattopadhyay, SINP, INDIA
  22 //
  23 //Email :  indra.das@saha.ac.in
  24 //         sukalyan.chattopadhyay@saha.ac.in
  25 //
  26 // This class implements a hit reconstruction algorithm for the dimuon
  27 // high level trigger.
  28 // The algorithm finds 3 pad clusters by looking for unique pads with a charge
  29 // above a certain threshold. A centre of gravity type calculation is applied
  30 // to the three pads forming the cluster to find the hit's X or Y coordinate
  31 // along the non-bending and bending planes individually.
  32 // The sepperate X and Y coordinates are then merged to give the full coordinate
  33 // of the hit.
  34 /////////////////////////////////////////////////
  35
  36 #include "AliHLTMUONHitReconstructor.h"
  37 #include "AliHLTMUONRecHitsBlockStruct.h"
  38 #include <cstring>
  39
  40
  41 const int AliHLTMUONHitReconstructor::fgkDetectorId = 0xA00;
  42 const int AliHLTMUONHitReconstructor::fgkDDLOffSet = 12 ;
  43 const int AliHLTMUONHitReconstructor::fgkNofDDL = 8 ;
  44
  45 const int AliHLTMUONHitReconstructor::fgkDDLHeaderSize = 8;
  46
  47 const int AliHLTMUONHitReconstructor::fgkEvenLutSize = 1645632 + 1;
  48 const int AliHLTMUONHitReconstructor::fgkOddLutSize = 3363840 + 1;
  49
  50 const int AliHLTMUONHitReconstructor::fgkLutLine[2] = {54208, 59648};
  51
  52 const int AliHLTMUONHitReconstructor::fgkMinIdManuChannel[2] = {917696, 64};
  53 const int AliHLTMUONHitReconstructor::fgkMaxIdManuChannel[2] = {2563327, 3363903};
  54
  55 const float AliHLTMUONHitReconstructor::fgkHalfPadSize[3] = {1.25, 2.50, 5.00};
  56
  57
  58 AliHLTMUONHitReconstructor::AliHLTMUONHitReconstructor():
  59   fkBlockHeaderSize(8),
  60   fkDspHeaderSize(8),
  61   fkBuspatchHeaderSize(4),
  62   fDCCut(0),
  63   fPadData(NULL),
  64   fLookUpTableData(NULL),
  65   fRecPoints(NULL),
  66   fRecPointsCount(NULL),
  67   fMaxRecPointsCount(0),
  68   fCentralCountB(0),
  69   fCentralCountNB(0),
  70   fIdOffSet(0),
  71   fDDLId(0),
  72   fDigitPerDDL(0),
  73   fDetManuChannelIdList(NULL),
  74   fCentralChargeB(NULL),
  75   fCentralChargeNB(NULL),
  76   fRecX(NULL),
  77   fRecY(NULL),
  78   fAvgChargeX(NULL),
  79   fAvgChargeY(NULL),
  80   fNofFiredDetElem(0),
  81   fDebugLevel(0),
  82   fBusToDetElem(),
  83   fBusToDDL()
  84 {
  85   // ctor
  86
  87   if(AliHLTMUONHitReconstructor::fgkEvenLutSize > AliHLTMUONHitReconstructor::fgkOddLutSize){
  88     fPadData = new DHLTPad[AliHLTMUONHitReconstructor::fgkEvenLutSize];
  89   }
  90   else{
  91     fPadData = new DHLTPad[AliHLTMUONHitReconstructor::fgkOddLutSize];
  92   }
  93
  94   fkBlockHeaderSize    = 8;
  95   fkDspHeaderSize      = 8;
  96   fkBuspatchHeaderSize = 4;
  97
  98   bzero(fGetIdTotalData,336*80*2*sizeof(int));
  99 }
 100
 101
 102 AliHLTMUONHitReconstructor::~AliHLTMUONHitReconstructor()
 103 {
 104   // dtor
 105
 106   if(fPadData){
 107     delete []fPadData;
 108     fPadData = NULL;
 109   }
 110
 111   if(fLookUpTableData){
 112     delete []fLookUpTableData;
 113     fLookUpTableData = NULL;
 114   }
 115 }
 116
 117
 118 int AliHLTMUONHitReconstructor::GetLutLine(int iDDL) const
 119 {
 120         return ( iDDL<16 ) ? fgkLutLine[0] : fgkLutLine[1] ;
 121 }
 122
 123
 124 bool AliHLTMUONHitReconstructor::LoadLookUpTable(DHLTLut* lookUpTableData, int lookUpTableId)
 125 {
 126   // function that loads LookUpTable (= position of each pad with electronic channel associated with it)
 127
 128   if(lookUpTableId<fgkDDLOffSet || lookUpTableId>= fgkDDLOffSet + fgkNofDDL){
 129     HLTError("DDL number is out of range (must be %d<=iDDL<%d)\n",fgkDDLOffSet,fgkDDLOffSet+fgkNofDDL);
 130     return false;
 131   }
 132
 133   fDDLId = lookUpTableId;
 134
 135   int lutSize = ((lookUpTableId%2)==0) ? fgkEvenLutSize : fgkOddLutSize ;
 136   int nofLutLine = GetLutLine(lookUpTableId);
 137   int idOffSet = fgkMinIdManuChannel[lookUpTableId%2];
 138
 139   int detManuChannelId;
 140
 141   fLookUpTableData = new DHLTLut[lutSize];
 142
 143   fLookUpTableData[0].fIdManuChannel = 0;
 144   fLookUpTableData[0].fIX = 0 ;
 145   fLookUpTableData[0].fIY = 0 ;
 146   fLookUpTableData[0].fRealX = 0.0 ;
 147   fLookUpTableData[0].fRealY = 0.0 ;
 148   fLookUpTableData[0].fRealZ = 0.0 ;
 149   fLookUpTableData[0].fPlane = -1 ;
 150   fLookUpTableData[0].fPcbZone = -1 ;
 151
 152   for(int i=0; i<nofLutLine; i++){
 153
 154     detManuChannelId = lookUpTableData[i].fIdManuChannel - idOffSet + 1;
 155     fLookUpTableData[detManuChannelId].fIdManuChannel = lookUpTableData[i].fIdManuChannel - idOffSet;
 156     fLookUpTableData[detManuChannelId].fIX = lookUpTableData[i].fIX ;
 157     fLookUpTableData[detManuChannelId].fIY = lookUpTableData[i].fIY ;
 158     fLookUpTableData[detManuChannelId].fRealX = lookUpTableData[i].fRealX ;
 159     fLookUpTableData[detManuChannelId].fRealY = lookUpTableData[i].fRealY ;
 160     fLookUpTableData[detManuChannelId].fRealZ = lookUpTableData[i].fRealZ ;
 161     fLookUpTableData[detManuChannelId].fPcbZone = lookUpTableData[i].fPcbZone ;
 162     fLookUpTableData[detManuChannelId].fPlane = lookUpTableData[i].fPlane ;
 163   }
 164   return true;
 165 }
 166
 167
 168 bool AliHLTMUONHitReconstructor::SetBusToDetMap(BusToDetElem busToDetElem)
 169 {
 170
 171   // function that loads BusPatch To Detection Element (SlatId) map
 172
 173   if(busToDetElem.size()==0) {
 174     HLTError("Empty BusToDetElem mapping");
 175     return false;
 176   } else {
 177     fBusToDetElem = busToDetElem;
 178   }
 179
 180   return true;
 181 }
 182
 183
 184 bool AliHLTMUONHitReconstructor::SetBusToDDLMap(BusToDDL busToDDL)
 185 {
 186
 187   // function that loads BusPatch To DDL Element (DDL) map
 188
 189   if(busToDDL.size()==0) {
 190     HLTError("Empty BusToDDL mapping");
 191     return false;
 192   } else {
 193     fBusToDDL = busToDDL;
 194   }
 195
 196   return true;
 197 }
 198
 199
 200 bool AliHLTMUONHitReconstructor::Run(
 201                 const AliHLTUInt32_t* rawData,
 202                 AliHLTUInt32_t rawDataSize,
 203                 AliHLTMUONRecHitStruct* recHit,
 204                 AliHLTUInt32_t& nofHit
 205         )
 206 {
 207   // main function called by HLTReconstructor to perform DHLT Hitreconstruction
 208
 209   fRecPoints = recHit;
 210   fMaxRecPointsCount = nofHit;
 211   fRecPointsCount = &nofHit;
 212   *fRecPointsCount = 0;
 213
 214   fPadData[0].fDetElemId = 0;
 215   fPadData[0].fBuspatchId = 0;
 216   fPadData[0].fIdManuChannel = 0;
 217   fPadData[0].fIX = 0 ;
 218   fPadData[0].fIY = 0 ;
 219   fPadData[0].fRealX = 0.0 ;
 220   fPadData[0].fRealY = 0.0 ;
 221   fPadData[0].fRealZ = 0.0 ;
 222   fPadData[0].fPlane = -1 ;
 223   fPadData[0].fPcbZone = -1 ;
 224   fPadData[0].fCharge = 0 ;
 225
 226   if (not ReadDDL(rawData, rawDataSize)) {
 227     // Dont need to log any message again. Already done so in ReadDDL.
 228     Clear();
 229     return false;
 230   }
 231
 232   if (fDigitPerDDL == 0)
 233   {
 234     // There are no digits to process so stop here.
 235     return true;
 236   }
 237
 238   if (not FindRecHits()) {
 239     HLTError("Failed to generate RecHits");
 240     Clear();
 241     return false;
 242   }
 243
 244   return true;
 245 }
 246
 247
 248 bool AliHLTMUONHitReconstructor::ReadDDL(
 249                 const AliHLTUInt32_t* rawData,
 250                 AliHLTUInt32_t rawDataSize
 251         )
 252 {
 253   //function to read Raw Data files
 254
 255   //const int* buffer = reinterpret_cast<const int*>(rawData);
 256   const AliHLTUInt32_t* buffer = rawData;
 257
 258   fIdOffSet= fgkMinIdManuChannel[(fDDLId%2)];
 259   fDetManuChannelIdList = new int[rawDataSize];
 260
 261   int index = 0;
 262   int dataCount = 0;
 263   fNofFiredDetElem = 0;
 264   int detElemId = 0 ;
 265   int buspatchId = 0;
 266   int prevDetElemId = 0 ;
 267   int totalBlockSize,blockRawDataSize;
 268   int totalDspSize,dspRawDataSize;
 269   int totalBuspatchSize,buspatchRawDataSize;
 270   int indexDsp,indexBuspatch,indexRawData;
 271   unsigned int dataWord;
 272   int charge;
 273   int idManuChannel;
 274
 275   for(int iBlock = 0; iBlock < 2 ;iBlock++){  // loop over 2 blocks
 276     totalBlockSize = buffer[index + 1];
 277     blockRawDataSize = buffer[index + 2];
 278     indexDsp = index + fkBlockHeaderSize;
 279     HLTDebug("totalBlockSize : %d, blockRawDataSize : %d",totalBlockSize,blockRawDataSize);
 280     while(blockRawDataSize > 0){
 281       totalDspSize = buffer[indexDsp + 1];
 282       dspRawDataSize = buffer[indexDsp + 2];
 283       HLTDebug("   totalDSPSize : %d, dspRawDataSize : %d",totalDspSize,dspRawDataSize);
 284       //if(buffer[indexDsp+1] == 1)
 285       dspRawDataSize --;                              // temporary solution to read buspatches
 286       indexBuspatch = indexDsp + fkDspHeaderSize + 2; // this extra 2 word comes from the faulty defination of Dsp header size
 287       while(dspRawDataSize > 0){
 288         totalBuspatchSize = buffer[indexBuspatch + 1];
 289         buspatchRawDataSize = buffer[indexBuspatch + 2];
 290         buspatchId = buffer[indexBuspatch + 3];
 291         if((detElemId = fBusToDetElem[buspatchId])==0){
 292           HLTError("No Detection element found for buspatch : %d",buspatchId);
 293           return false;
 294         }
 295         HLTDebug("\ttotalBusPatchSize : %d, buspatchRawDataSize : %d",totalBuspatchSize,buspatchRawDataSize);
 296         indexRawData = indexBuspatch + fkBuspatchHeaderSize;
 297         while(buspatchRawDataSize > 0){
 298           dataWord = buffer[indexRawData];
 299           charge = (unsigned short)(dataWord & 0xFFF);
 300
 301           idManuChannel = 0x0;
 302           idManuChannel = (idManuChannel|(detElemId%100))<<17;
 303           idManuChannel |= (dataWord >> 12) & 0x1FFFF;
 304           idManuChannel -= fIdOffSet ;
 305
 306           if(charge > fDCCut && charge > 5){  // (charge > 4) is due cut out the noise level
 307             fPadData[idManuChannel].fBuspatchId = buspatchId;
 308             fPadData[idManuChannel].fDetElemId = detElemId;
 309             fPadData[idManuChannel].fIdManuChannel = idManuChannel;
 310             fPadData[idManuChannel].fIX = fLookUpTableData[idManuChannel+1].fIX;
 311             fPadData[idManuChannel].fIY = fLookUpTableData[idManuChannel+1].fIY;
 312             fPadData[idManuChannel].fRealX = fLookUpTableData[idManuChannel+1].fRealX;
 313             fPadData[idManuChannel].fRealY = fLookUpTableData[idManuChannel+1].fRealY;
 314             fPadData[idManuChannel].fRealZ = fLookUpTableData[idManuChannel+1].fRealZ;
 315             fPadData[idManuChannel].fPcbZone = fLookUpTableData[idManuChannel+1].fPcbZone;
 316             fPadData[idManuChannel].fPlane = fLookUpTableData[idManuChannel+1].fPlane;
 317             fPadData[idManuChannel].fCharge = charge;
 318
 319             fDetManuChannelIdList[dataCount] = idManuChannel;
 320             if(detElemId != prevDetElemId){
 321               if(fNofFiredDetElem>0){
 322                 fMaxFiredPerDetElem[fNofFiredDetElem-1] = dataCount;
 323               }
 324               fNofFiredDetElem++;
 325               prevDetElemId = detElemId ;
 326             }
 327
 328             HLTDebug("\t  buspatch : %d, detele : %d, id : %d, manu : %d, channel : %d, X : %f, Y: %f",
 329                     fPadData[idManuChannel].fBuspatchId,fPadData[idManuChannel].fDetElemId,
 330                     idManuChannel,((dataWord >> 12) & 0x7FF),((dataWord >> 23) & 0x3F),
 331                     fPadData[idManuChannel].fRealX,fPadData[idManuChannel].fRealY);
 332
 333             dataCount ++;
 334           }// if charge is more than DC Cut limit condition
 335
 336           indexRawData++;
 337           buspatchRawDataSize --;
 338         }
 339         indexBuspatch += totalBuspatchSize;
 340         dspRawDataSize -= totalBuspatchSize;
 341       }// buspatch loop
 342       indexDsp += totalDspSize;
 343       blockRawDataSize -= totalDspSize;
 344     }// DSP loop
 345     index = totalBlockSize;
 346   }// Block loop
 347
 348   fDigitPerDDL = dataCount;
 349   fMaxFiredPerDetElem[fNofFiredDetElem-1] = dataCount;
 350
 351   if(fDigitPerDDL == 0){
 352     HLTInfo("An Empty DDL File found");
 353   }
 354
 355   return true;
 356 }
 357
 358
 359 bool AliHLTMUONHitReconstructor::FindRecHits()
 360 {
 361   // fuction that calls hit reconstruction detector element-wise
 362
 363   for(int iDet=0; iDet<fNofFiredDetElem ; iDet++){
 364
 365     fCentralCountB = 0 ;
 366     fCentralCountNB = 0 ;
 367     fCentralChargeB = new int[fMaxFiredPerDetElem[iDet]];
 368     fCentralChargeNB = new int[fMaxFiredPerDetElem[iDet]];
 369
 370     if(iDet>0)
 371       FindCentralHits(fMaxFiredPerDetElem[iDet-1],fMaxFiredPerDetElem[iDet]);
 372     else
 373       FindCentralHits(0,fMaxFiredPerDetElem[iDet]);
 374
 375     RecXRecY();
 376     if(!MergeRecHits()){
 377       HLTError("Failed to merge hits\n");
 378       return false;
 379     }
 380
 381     if(iDet==0)
 382       for(int i=0;i<fMaxFiredPerDetElem[iDet];i++)
 383         fGetIdTotalData[fPadData[fDetManuChannelIdList[i]].fIX][fPadData[fDetManuChannelIdList[i]].fIY][fPadData[fDetManuChannelIdList[i]].fPlane] = 0;
 384     else
 385       for(int i=fMaxFiredPerDetElem[iDet-1];i<fMaxFiredPerDetElem[iDet];i++)
 386         fGetIdTotalData[fPadData[fDetManuChannelIdList[i]].fIX][fPadData[fDetManuChannelIdList[i]].fIY][fPadData[fDetManuChannelIdList[i]].fPlane] = 0;
 387
 388     if(fCentralChargeB){
 389       delete []fCentralChargeB;
 390       fCentralChargeB = NULL;
 391     }
 392
 393     if(fCentralChargeNB){
 394       delete []fCentralChargeNB;
 395       fCentralChargeNB = NULL;
 396     }
 397
 398   }
 399
 400   Clear();
 401
 402   return true;
 403 }
 404
 405
 406 void AliHLTMUONHitReconstructor::FindCentralHits(int minPadId, int maxPadId)
 407 {
 408   // to find central hit associated with each cluster
 409
 410   int b,nb;
 411   int idManuChannelCentral;
 412   bool hasFind;
 413   int idManuChannel;
 414
 415   for(int iPad=minPadId;iPad<maxPadId;iPad++){
 416     idManuChannel   = fDetManuChannelIdList[iPad];
 417
 418
 419     fGetIdTotalData[fPadData[idManuChannel].fIX]
 420       [fPadData[idManuChannel].fIY]
 421       [fPadData[idManuChannel].fPlane] = idManuChannel;
 422
 423     if(fPadData[idManuChannel].fPlane == 0 ){//&& fPadData[idManuChannel].fIY > (0+1) && fPadData[idManuChannel].fIY < (79 - 1)){
 424       //if(fPadData[idManuChannel].fIY > 0){
 425       if(fCentralCountB>0){
 426         hasFind = false;
 427         for(b = 0;b<fCentralCountB;b++){
 428           idManuChannelCentral = fCentralChargeB[b];
 429           if(fPadData[idManuChannel].fIX == fPadData[idManuChannelCentral].fIX
 430              &&
 431              (fPadData[idManuChannel].fIY
 432               == fPadData[idManuChannelCentral].fIY + 1
 433               ||
 434               fPadData[idManuChannel].fIY
 435               == fPadData[idManuChannelCentral].fIY + 2
 436               ||
 437               fPadData[idManuChannel].fIY
 438               == fPadData[idManuChannelCentral].fIY - 2
 439               ||
 440               fPadData[idManuChannel].fIY
 441               == fPadData[idManuChannelCentral].fIY - 1)){
 442
 443             hasFind = true;
 444             if(fPadData[idManuChannel].fCharge > fPadData[idManuChannelCentral].fCharge){
 445               fCentralChargeB[b] = idManuChannel;
 446             }// if condn on pad charge
 447           }// if condon on pad position
 448         }// for loop over b
 449         if(!hasFind){
 450           fCentralChargeB[fCentralCountB] = idManuChannel;
 451           fCentralCountB++;
 452         }
 453       }
 454       else{
 455         fCentralChargeB[fCentralCountB] = idManuChannel;
 456         fCentralCountB++;
 457       }// check the size of centralHitB
 458       for(b = 0;b<fCentralCountB;b++){
 459         idManuChannelCentral = fCentralChargeB[b];
 460       }
 461       //}// if cond on iY > 2 (to avoid edge value pb)
 462     }// B/Nb checking
 463     else{
 464       if(fCentralCountNB>0){
 465         hasFind = false;
 466         for(nb = 0;nb<fCentralCountNB;nb++){
 467           idManuChannelCentral = fCentralChargeNB[nb];
 468           if(fPadData[idManuChannel].fIY == fPadData[idManuChannelCentral].fIY
 469              &&
 470              (fPadData[idManuChannel].fIX
 471               == fPadData[idManuChannelCentral].fIX + 1
 472               ||
 473               fPadData[idManuChannel].fIX
 474               == fPadData[idManuChannelCentral].fIX + 2
 475               ||
 476               fPadData[idManuChannel].fIX
 477               == fPadData[idManuChannelCentral].fIX - 2
 478               ||
 479               fPadData[idManuChannel].fIX
 480               == fPadData[idManuChannelCentral].fIX - 1)){
 481
 482             hasFind = true;
 483             if(fPadData[idManuChannel].fCharge > fPadData[idManuChannelCentral].fCharge){
 484               fCentralChargeNB[nb] = idManuChannel;
 485             }// if condn over to find higher charge
 486           }// if condn over to find position
 487         }// for loop over presently all nb values
 488         if(!hasFind){
 489           fCentralChargeNB[fCentralCountNB] = idManuChannel;
 490           fCentralCountNB++;
 491         }
 492       }// centralHitNB size test
 493       else{
 494         fCentralChargeNB[fCentralCountNB] = idManuChannel;
 495         fCentralCountNB++;
 496       }// centralHitNB size test
 497
 498     }// fill for bending and nonbending hit
 499   }// detElemId loop
 500 }
 501
 502
 503 void AliHLTMUONHitReconstructor::RecXRecY()
 504 {
 505   // find reconstructed X and Y for each plane separately
 506   int b,nb;
 507   int idCentral;
 508   int idLower = 0;
 509   int idUpper = 0;
 510   int idRight = 0;
 511   int idLeft = 0;
 512   fRecY = new float[fCentralCountB];
 513   fRecX = new float[fCentralCountNB];
 514
 515   fAvgChargeY = new float[fCentralCountB];
 516   fAvgChargeX = new float[fCentralCountNB];
 517
 518   for(b=0;b<fCentralCountB;b++){
 519     idCentral = fCentralChargeB[b];
 520
 521     if(fPadData[idCentral].fIY==0)
 522       idLower = 0;
 523     else
 524       idLower = fGetIdTotalData[fPadData[idCentral].fIX][fPadData[idCentral].fIY-1][0];
 525
 526     if(fPadData[idCentral].fIY==79)
 527       idUpper = 0;
 528     else
 529       idUpper = fGetIdTotalData[fPadData[idCentral].fIX][fPadData[idCentral].fIY+1][0];
 530
 531     fRecY[b] = (fPadData[idCentral].fRealY*fPadData[idCentral].fCharge
 532                +
 533                fPadData[idUpper].fRealY*fPadData[idUpper].fCharge
 534                +
 535                fPadData[idLower].fRealY*fPadData[idLower].fCharge
 536                )/(fPadData[idCentral].fCharge + fPadData[idUpper].fCharge + fPadData[idLower].fCharge) ;
 537
 538     fAvgChargeY[b] = fPadData[idCentral].fCharge;
 539
 540
 541     fRecY[b] += 0.025*sin(12.56637*(0.25-(fRecY[b] - fPadData[idCentral].fRealY))) ;
 542   }
 543
 544   for(nb=0;nb<fCentralCountNB;nb++){
 545     idCentral = fCentralChargeNB[nb];
 546
 547     if(fPadData[idCentral].fIX==0)
 548       idLeft = 0;
 549     else
 550       idLeft = fGetIdTotalData[fPadData[idCentral].fIX-1][fPadData[idCentral].fIY][1];
 551
 552     if(fPadData[idCentral].fIX==335)
 553       idRight = 0 ;
 554     else
 555       idRight = fGetIdTotalData[fPadData[idCentral].fIX+1][fPadData[idCentral].fIY][1];
 556
 557     fRecX[nb] = (fPadData[idCentral].fRealX*fPadData[idCentral].fCharge
 558                  +
 559                  fPadData[idRight].fRealX*fPadData[idRight].fCharge
 560                  +
 561                  fPadData[idLeft].fRealX*fPadData[idLeft].fCharge
 562                  )/(fPadData[idCentral].fCharge + fPadData[idRight].fCharge + fPadData[idLeft].fCharge);
 563
 564
 565     fAvgChargeX[nb] = fPadData[idCentral].fCharge;
 566
 567   }
 568 }
 569
 570
 571 bool AliHLTMUONHitReconstructor::MergeRecHits()
 572 {
 573   // Merge reconstructed hits first over same plane then bending plane with non-bending plane
 574
 575   int idCentralB,idCentralNB ;
 576   float padCenterXB;
 577   float padCenterYNB;
 578   float diffX,diffY;
 579   float halfPadLengthX,halfPadLengthY;
 580
 581   // MERGE Bending Plane hits, which are placed side by side
 582   for(int i=0;i<fCentralCountB-1;i++){
 583     if(fRecY[i] != 0.0){
 584       for(int j=i+1;j<fCentralCountB;j++){
 585
 586         if(fCentralChargeB[i]==fCentralChargeB[j]){
 587           fRecY[j] = 0.0;
 588           continue;
 589         }
 590         else if(
 591            (
 592             fPadData[fCentralChargeB[i]].fIY == fPadData[fCentralChargeB[j]].fIY
 593             )
 594            &&
 595            (
 596             fPadData[fCentralChargeB[i]].fIX == fPadData[fCentralChargeB[j]].fIX + 1
 597             ||
 598             fPadData[fCentralChargeB[i]].fIX == fPadData[fCentralChargeB[j]].fIX - 1
 599             )
 600            &&
 601            fRecY[j] != 0.0
 602            &&
 603            fRecY[i] != 0.0
 604            ){
 605
 606           if(fAvgChargeY[i] > fAvgChargeY[j]){
 607             fRecY[i] = (fRecY[i]*fAvgChargeY[i] + fRecY[j]*fAvgChargeY[j]
 608                         )/(fAvgChargeY[i] + fAvgChargeY[j]);
 609             fRecY[j] = 0.0;
 610           }
 611           else{
 612             fRecY[j] = (fRecY[i]*fAvgChargeY[i] + fRecY[j]*fAvgChargeY[j]
 613                         )/(fAvgChargeY[i] + fAvgChargeY[j]);
 614             fRecY[i] = 0.0;
 615
 616           }// search for higher charge
 617         }//pad position
 618       }//j for loop
 619     }//if fRecY[i] != 0.0
 620   }// i for loop
 621
 622   // MERGE Non Bending Plane hits, which are place side by side
 623   for(int i=0;i<fCentralCountNB-1;i++){
 624     if(fRecX[i] != 0.0){
 625       for(int j=i+1;j<fCentralCountNB;j++){
 626
 627         if(fCentralChargeNB[i]==fCentralChargeNB[j]){
 628           fRecX[j] = 0.0;
 629           continue;
 630         }
 631         else if(
 632            (
 633             fPadData[fCentralChargeNB[i]].fIX == fPadData[fCentralChargeNB[j]].fIX
 634             )
 635            &&
 636            (
 637             fPadData[fCentralChargeNB[i]].fIY == fPadData[fCentralChargeNB[j]].fIY + 1
 638             ||
 639             fPadData[fCentralChargeNB[i]].fIY == fPadData[fCentralChargeNB[j]].fIY - 1
 640             )
 641            &&
 642            fRecX[j] != 0.0
 643            &&
 644            fRecX[i] != 0.0
 645            ){
 646
 647           if(fAvgChargeX[i] > fAvgChargeX[j]){
 648             fRecX[i] = (fRecX[i]*fAvgChargeX[i] + fRecX[j]*fAvgChargeX[j]
 649                        )/(fAvgChargeX[i] + fAvgChargeX[j]);
 650             fRecX[j] = 0.0;
 651           }
 652           else{
 653             fRecX[j] = (fRecX[i]*fAvgChargeX[i] + fRecX[j]*fAvgChargeX[j]
 654                        )/(fAvgChargeX[i] + fAvgChargeX[j]);
 655             fRecX[i] = 0.0;
 656           }// search for higher charge
 657         }//pad position
 658       }//j for loop
 659     }//if fRecX[i] != 0.0
 660   }// i for loop
 661
 662   // Merge bending Plane hits with Non Bending
 663   for(int b=0;b<fCentralCountB;b++){
 664     if(fRecY[b]!=0.0){
 665       idCentralB = fCentralChargeB[b];
 666       padCenterXB = fPadData[idCentralB].fRealX;
 667
 668       halfPadLengthX = fgkHalfPadSize[fPadData[idCentralB].fPcbZone] ;
 669
 670       for(int nb=0;nb<fCentralCountNB;nb++){
 671         if(fRecX[nb]!=0.0){
 672           idCentralNB = fCentralChargeNB[nb];
 673
 674           padCenterYNB = fPadData[idCentralNB].fRealY;
 675
 676           halfPadLengthY = fgkHalfPadSize[fPadData[idCentralNB].fPcbZone] ;
 677
 678           if(fabsf(fRecX[nb]) > fabsf(padCenterXB))
 679             diffX = fabsf(fRecX[nb]) -  fabsf(padCenterXB);
 680           else
 681             diffX = fabsf(padCenterXB) -  fabsf(fRecX[nb]);
 682
 683           if(fabsf(padCenterYNB)>fabsf(fRecY[b]))
 684             diffY = fabsf(padCenterYNB) - fabsf(fRecY[b]);
 685           else
 686             diffY =  fabsf(fRecY[b]) - fabsf(padCenterYNB);
 687
 688           if(diffX < halfPadLengthX && diffY < halfPadLengthY ){//&& fPadData[idCentralB].fIY != 0){
 689
 690             // First check that we have not overflowed the buffer.
 691             if((*fRecPointsCount) == fMaxRecPointsCount){
 692               HLTError("Nof RecHit (i.e. %d) exceeds the max nof RecHit limit %d\n",(*fRecPointsCount),fMaxRecPointsCount);
 693               return false;
 694             }
 695
 696             //fRecPoints[(*fRecPointsCount)].fId = idCentralB;
 697             fRecPoints[(*fRecPointsCount)].fX = fRecX[nb];
 698             fRecPoints[(*fRecPointsCount)].fY = fRecY[b];
 699             fRecPoints[(*fRecPointsCount)].fZ = fPadData[idCentralB].fRealZ;
 700             //fRecPoints[(*fRecPointsCount)].fDetElemId = (AliHLTUInt32_t)fPadData[idCentralB].fDetElemId;
 701             (*fRecPointsCount)++;
 702           }//if lies wihtin 5.0 mm
 703         }// condn over fRecX ! = 0.0
 704       }// loop over NB side
 705     }// condn on fRecY[b] !=  0.0
 706   }// loop over B side;
 707
 708   if(fRecX){
 709     delete []fRecX;
 710     fRecX = NULL;
 711   }
 712
 713   if(fRecY){
 714     delete []fRecY;
 715     fRecY = NULL;
 716   }
 717
 718   if(fAvgChargeX){
 719     delete []fAvgChargeX;
 720     fAvgChargeX = NULL;
 721   }
 722
 723   if(fAvgChargeY){
 724     delete []fAvgChargeY;
 725     fAvgChargeY = NULL;
 726   }
 727
 728   return true;
 729 }
 730
 731
 732 void AliHLTMUONHitReconstructor::Clear()
 733 {
 734   for(int iPad=0;iPad<fDigitPerDDL;iPad++){
 735     fGetIdTotalData[fPadData[fDetManuChannelIdList[iPad]].fIX][fPadData[fDetManuChannelIdList[iPad]].fIY][fPadData[fDetManuChannelIdList[iPad]].fPlane] = 0;
 736     fPadData[fDetManuChannelIdList[iPad]].fDetElemId = 0;
 737     fPadData[fDetManuChannelIdList[iPad]].fBuspatchId = 0;
 738     fPadData[fDetManuChannelIdList[iPad]].fIdManuChannel = 0;
 739     fPadData[fDetManuChannelIdList[iPad]].fIX = 0 ;
 740     fPadData[fDetManuChannelIdList[iPad]].fIY = 0 ;
 741     fPadData[fDetManuChannelIdList[iPad]].fRealX = 0.0 ;
 742     fPadData[fDetManuChannelIdList[iPad]].fRealY = 0.0 ;
 743     fPadData[fDetManuChannelIdList[iPad]].fRealZ = 0.0 ;
 744     fPadData[fDetManuChannelIdList[iPad]].fPlane = -1 ;
 745     fPadData[fDetManuChannelIdList[iPad]].fPcbZone = -1 ;
 746     fPadData[fDetManuChannelIdList[iPad]].fCharge = 0 ;
 747   }
 748
 749   for(int i=0;i<13;i++)
 750     fMaxFiredPerDetElem[i] = 0;
 751
 752   if(fDetManuChannelIdList){
 753     delete []fDetManuChannelIdList;
 754     fDetManuChannelIdList = NULL;
 755   }
 756
 757   if(fCentralChargeB){
 758     delete []fCentralChargeB;
 759     fCentralChargeB = NULL;
 760   }
 761
 762   if(fCentralChargeNB){
 763     delete []fCentralChargeNB;
 764     fCentralChargeNB = NULL;
 765   }
 766
 767   if(fRecX){
 768     delete []fRecX;
 769     fRecX = NULL;
 770   }
 771
 772   if(fRecY){
 773     delete []fRecY;
 774     fRecY = NULL;
 775   }
 776
 777   if(fAvgChargeX){
 778     delete []fAvgChargeX;
 779     fAvgChargeX = NULL;
 780   }
 781
 782   if(fAvgChargeY){
 783     delete []fAvgChargeY;
 784     fAvgChargeY = NULL;
 785   }
 786 }