HLT/hough/AliL3HoughMaxFinder.cxx

   1 #include <string.h>
   2
   3 #include <TH2.h>
   4
   5 #include "AliL3TrackArray.h"
   6 #include "AliL3HoughTrack.h"
   7 #include "AliL3HoughMaxFinder.h"
   8
   9 ClassImp(AliL3HoughMaxFinder)
  10
  11
  12 AliL3HoughMaxFinder::AliL3HoughMaxFinder()
  13 {
  14   //Default constructor
  15   fThreshold = 0;
  16   //fTracks = 0;
  17   fHistoType=0;
  18 }
  19
  20
  21 AliL3HoughMaxFinder::AliL3HoughMaxFinder(Char_t *histotype)
  22 {
  23   //Constructor
  24
  25   //fTracks = new AliL3TrackArray("AliL3HoughTrack");
  26   if(strcmp(histotype,"KappaPhi")==0) fHistoType='c';
  27   if(strcmp(histotype,"DPsi")==0) fHistoType='l';
  28
  29 }
  30
  31
  32 AliL3HoughMaxFinder::~AliL3HoughMaxFinder()
  33 {
  34   //Destructor
  35
  36 }
  37
  38 AliL3TrackArray *AliL3HoughMaxFinder::FindBigMaxima(TH2F *hist)
  39 {
  40   Int_t xmin = hist->GetXaxis()->GetFirst();
  41   Int_t xmax = hist->GetXaxis()->GetLast();
  42   Int_t ymin = hist->GetYaxis()->GetFirst();
  43   Int_t ymax = hist->GetYaxis()->GetLast();
  44   Int_t bin[25],bin_index;
  45   Stat_t value[25];
  46
  47   AliL3TrackArray *tracks = new AliL3TrackArray("AliL3HoughTrack");
  48   AliL3HoughTrack *track;
  49
  50   for(Int_t xbin=xmin+2; xbin<xmax-3; xbin++)
  51     {
  52       for(Int_t ybin=ymin+2; ybin<ymax-3; ybin++)
  53         {
  54           bin_index=0;
  55           for(Int_t xb=xbin-2; xb<xbin+3; xb++)
  56             {
  57               for(Int_t yb=ybin-2; yb<ybin+3; yb++)
  58                 {
  59                   bin[bin_index]=hist->GetBin(xb,yb);
  60                   value[bin_index]=hist->GetBinContent(bin[bin_index]);
  61                   bin_index++;
  62                 }
  63
  64             }
  65           if(value[12]==0) continue;
  66           Int_t b=0;
  67           while(1)
  68             {
  69               if(value[b] > value[12] || b==bin_index) break;
  70               b++;
  71               printf("b %d\n",b);
  72             }
  73           if(b == bin_index)
  74             {
  75               //Found maxima
  76               Double_t max_x = hist->GetXaxis()->GetBinCenter(xbin);
  77               Double_t max_y = hist->GetYaxis()->GetBinCenter(ybin);
  78               track = (AliL3HoughTrack*)tracks->NextTrack();
  79               track->SetTrackParameters(max_x,max_y,value[12]);
  80             }
  81         }
  82     }
  83
  84   tracks->QSort();
  85   return tracks;
  86 }
  87
  88
  89 AliL3TrackArray *AliL3HoughMaxFinder::FindMaxima(TH2F *hist,Int_t *rowrange,Int_t ref_row)
  90 {
  91   //Locate all the maxima in input histogram.
  92   //Maxima is defined as bins with more entries than the
  93   //immediately neighbouring bins.
  94
  95   Int_t xmin = hist->GetXaxis()->GetFirst();
  96   Int_t xmax = hist->GetXaxis()->GetLast();
  97   Int_t ymin = hist->GetYaxis()->GetFirst();
  98   Int_t ymax = hist->GetYaxis()->GetLast();
  99   Int_t bin[9],track_counter=0;
 100   Stat_t value[9];
 101
 102   AliL3TrackArray *tracks = new AliL3TrackArray("AliL3HoughTrack");
 103   AliL3HoughTrack *track;
 104
 105   for(Int_t xbin=xmin+1; xbin<xmax-1; xbin++)
 106     {
 107       for(Int_t ybin=ymin+1; ybin<ymax-1; ybin++)
 108         {
 109           bin[0] = hist->GetBin(xbin-1,ybin-1);
 110           bin[1] = hist->GetBin(xbin,ybin-1);
 111           bin[2] = hist->GetBin(xbin+1,ybin-1);
 112           bin[3] = hist->GetBin(xbin-1,ybin);
 113           bin[4] = hist->GetBin(xbin,ybin);
 114           bin[5] = hist->GetBin(xbin+1,ybin);
 115           bin[6] = hist->GetBin(xbin-1,ybin+1);
 116           bin[7] = hist->GetBin(xbin,ybin+1);
 117           bin[8] = hist->GetBin(xbin+1,ybin+1);
 118           value[0] = hist->GetBinContent(bin[0]);
 119           value[1] = hist->GetBinContent(bin[1]);
 120           value[2] = hist->GetBinContent(bin[2]);
 121           value[3] = hist->GetBinContent(bin[3]);
 122           value[4] = hist->GetBinContent(bin[4]);
 123           value[5] = hist->GetBinContent(bin[5]);
 124           value[6] = hist->GetBinContent(bin[6]);
 125           value[7] = hist->GetBinContent(bin[7]);
 126           value[8] = hist->GetBinContent(bin[8]);
 127
 128           if(value[4] <= fThreshold) continue;//central bin below threshold
 129
 130           if(value[4]>value[0] && value[4]>value[1] && value[4]>value[2]
 131              && value[4]>value[3] && value[4]>value[5] && value[4]>value[6]
 132              && value[4]>value[7] && value[4]>value[8])
 133             {
 134               //Found a local maxima
 135               Float_t max_x = hist->GetXaxis()->GetBinCenter(xbin);
 136               Float_t max_y = hist->GetYaxis()->GetBinCenter(ybin);
 137
 138               track = (AliL3HoughTrack*)tracks->NextTrack();
 139
 140
 141               switch(fHistoType)
 142                 {
 143                 case 'c':
 144                   track->SetTrackParameters(max_x,max_y,(Int_t)value[4]);
 145                   break;
 146                 case 'l':
 147                   track->SetLineParameters(max_x,max_y,(Int_t)value[4],rowrange,ref_row);
 148                   break;
 149                 default:
 150                   printf("AliL3HoughMaxFinder: Error in tracktype\n");
 151                 }
 152
 153               track_counter++;
 154
 155
 156             }
 157           else
 158             continue; //not a maxima
 159         }
 160     }
 161   tracks->QSort();
 162   return tracks;
 163 }
 164
 165
 166 AliL3TrackArray *AliL3HoughMaxFinder::LookForPeaks(TH2F *hist,Int_t nbins)
 167 {
 168
 169   AliL3TrackArray *tracks = new AliL3TrackArray("AliL3HoughTrack");
 170   AliL3HoughTrack *track;
 171   Int_t xmin = hist->GetXaxis()->GetFirst();
 172   Int_t xmax = hist->GetXaxis()->GetLast();
 173   Int_t ymin = hist->GetYaxis()->GetFirst();
 174   Int_t ymax = hist->GetYaxis()->GetLast();
 175
 176   Int_t weight_loc;
 177   Int_t candidate=0;
 178   for(Int_t xbin=xmin+nbins; xbin <= xmax - nbins; xbin++)
 179     {
 180       for(Int_t ybin=ymin+nbins; ybin <= ymax - nbins; ybin++)
 181         {
 182           weight_loc=0;
 183           for(Int_t xbin_loc = xbin-nbins; xbin_loc <= xbin+nbins; xbin_loc++)
 184             {
 185               for(Int_t ybin_loc = ybin-nbins; ybin_loc <= ybin+nbins; ybin_loc++)
 186                 {
 187                   Int_t bin_loc = hist->GetBin(xbin_loc,ybin_loc);
 188                   weight_loc += (Int_t)hist->GetBinContent(bin_loc);
 189                 }
 190             }
 191
 192           if(weight_loc > 0)
 193             {
 194               track = (AliL3HoughTrack*)tracks->NextTrack();
 195               track->SetTrackParameters(hist->GetXaxis()->GetBinCenter(xbin),hist->GetYaxis()->GetBinCenter(ybin),weight_loc);
 196             }
 197
 198         }
 199     }
 200   tracks->QSort();
 201
 202   AliL3HoughTrack *track1,*track2;
 203
 204   for(Int_t i=1; i<tracks->GetNTracks(); i++)
 205     {
 206       track1 = (AliL3HoughTrack*)tracks->GetCheckedTrack(i);
 207       if(!track1) continue;
 208       Int_t xbin1 = hist->GetXaxis()->FindBin(track1->GetKappa());
 209       Int_t ybin1 = hist->GetYaxis()->FindBin(track1->GetPhi0());
 210       for(Int_t j=0; j < i; j++)
 211         {
 212           track2 = (AliL3HoughTrack*)tracks->GetCheckedTrack(j);
 213           if(!track2) continue;
 214           Int_t xbin2 = hist->GetXaxis()->FindBin(track2->GetKappa());
 215           Int_t ybin2 = hist->GetYaxis()->FindBin(track2->GetPhi0());
 216           if(abs(xbin1-xbin2) < 10 && abs(ybin1-ybin2) < 10)
 217             {
 218               tracks->Remove(i);
 219               break;
 220             }
 221         }
 222
 223     }
 224   tracks->Compress();
 225   return tracks;
 226 }
 227
 228 AliL3TrackArray *AliL3HoughMaxFinder::LookInWindows(TH2F *hist,Int_t nbins,Int_t t1,Double_t t2,Int_t t3)
 229 {
 230
 231   AliL3TrackArray *tracks = new AliL3TrackArray("AliL3HoughTrack");
 232   AliL3HoughTrack *track = (AliL3HoughTrack*)tracks->NextTrack();
 233
 234   Int_t xmin = hist->GetXaxis()->GetFirst();
 235   Int_t xmax = hist->GetXaxis()->GetLast();
 236   Int_t ymin = hist->GetYaxis()->GetFirst();
 237   Int_t ymax = hist->GetYaxis()->GetLast();
 238
 239   for(Int_t xbin=xmin+nbins; xbin < xmax - nbins; xbin += nbins)
 240     {
 241       for(Int_t ybin=ymin+nbins; ybin<ymax-nbins; ybin += nbins)
 242         {
 243           //Int_t bin = hist->GetBin(xbin,ybin);
 244           //if((Int_t)hist->GetBinContent(bin)==0) continue;
 245
 246           Int_t xrange[2] = {xbin-nbins,xbin+nbins};
 247           Int_t yrange[2] = {ybin-nbins,ybin+nbins};
 248           if(LocatePeak(hist,track,xrange,yrange,t1,t2,t3))
 249             track = (AliL3HoughTrack*)tracks->NextTrack();
 250           else
 251             continue;
 252         }
 253     }
 254   tracks->QSort();
 255   tracks->RemoveLast();
 256   return tracks;
 257
 258 }
 259
 260 Bool_t AliL3HoughMaxFinder::LocatePeak(TH2F *hist,AliL3HoughTrack *track,Int_t *xrange,Int_t *yrange,Int_t t1,Double_t t2,Int_t t3)
 261 {
 262
 263   /*  Int_t xmin = hist->GetXaxis()->FindBin(track->GetKappa()) - nbins;
 264   Int_t xmax = hist->GetXaxis()->FindBin(track->GetKappa()) + nbins;
 265   Int_t ymin = hist->GetYaxis()->FindBin(track->GetPhi0()) - nbins;
 266   Int_t ymax = hist->GetYaxis()->FindBin(track->GetPhi0()) + nbins;
 267   Int_t nbinsx = nbins*2 + 1;
 268
 269   if(xmin < 0 || xmax > hist->GetXaxis()->GetNbins() || ymin < 0 || ymax > hist->GetYaxis()->GetNbins())
 270     return false;
 271   */
 272   Int_t xmin = xrange[0];
 273   Int_t xmax = xrange[1];
 274   Int_t ymin = yrange[0];
 275   Int_t ymax = yrange[1];
 276
 277   if(xmin < 0)
 278     xmin=0;
 279   if(xmax > hist->GetXaxis()->GetNbins())
 280     xmax = hist->GetXaxis()->GetNbins();
 281   if(ymin < 0)
 282     ymin=0;
 283   if(ymax > hist->GetYaxis()->GetNbins())
 284     ymax = hist->GetYaxis()->GetNbins();
 285
 286   printf("Defined xrange %d %d yrange %d %d\n",xmin,xmax,ymin,ymax);
 287
 288
 289   Int_t totbins = hist->GetXaxis()->GetNbins();
 290   Int_t *m = new Int_t[totbins];
 291   Int_t *m_low = new Int_t[totbins];
 292   Int_t *m_up = new Int_t[totbins];
 293
 294
 295   for(Int_t i=0; i<totbins; i++)
 296     {
 297       m[i]=0;
 298       m_low[i]=0;
 299       m_up[i]=0;
 300     }
 301
 302   Int_t max_x=0,sum=0,max_xbin=0,bin;
 303
 304   for(Int_t xbin=xmin; xbin<=xmax; xbin++)
 305     {
 306       for(Int_t ybin=ymin; ybin < ymax - t1; ybin++)
 307         {
 308           sum = 0;
 309           for(Int_t y=ybin; y < ybin+t1; y++)
 310             {
 311               //Inside window
 312               bin = hist->GetBin(xbin,y);
 313               sum += (Int_t)hist->GetBinContent(bin);
 314
 315             }
 316           if(sum > m[xbin]) //Max value locally in this xbin
 317             {
 318               m[xbin]=sum;
 319               m_low[xbin]=ybin;
 320               m_up[xbin]=ybin + t1 - 1;
 321             }
 322
 323         }
 324
 325       if(m[xbin] > max_x) //Max value globally in x-direction
 326         {
 327           max_xbin = xbin;
 328           max_x = m[xbin];//sum;
 329         }
 330     }
 331
 332   if(max_x == 0)
 333     {
 334       //printf("\nmax_x == 0\n");
 335       return false;
 336     }
 337   //printf("max_xbin %d max_x %d m_low %d m_up %d\n",max_xbin,max_x,m_low[max_xbin],m_up[max_xbin]);
 338   //printf("ylow %f yup %f\n",hist->GetYaxis()->GetBinCenter(m_low[max_xbin]),hist->GetYaxis()->GetBinCenter(m_up[max_xbin]));
 339
 340   //Determine a width in the x-direction
 341   Int_t x_low=0,x_up=0;
 342   for(Int_t xbin=max_xbin-1; xbin >= xmin; xbin--)
 343     {
 344       if(m[xbin]*t2 < max_x)
 345         {
 346           x_low = xbin+1;
 347           break;
 348         }
 349     }
 350   for(Int_t xbin = max_xbin+1; xbin <=xmax; xbin++)
 351     {
 352       if(m[xbin]*t2 < max_x)
 353         {
 354           x_up = xbin-1;
 355           break;
 356         }
 357     }
 358
 359   Double_t top=0,butt=0,value,x_peak;
 360
 361   // printf("xlow %f xup %f\n",hist->GetXaxis()->GetBinCenter(x_low),hist->GetXaxis()->GetBinCenter(x_up));
 362   //printf("Spread in x %d\n",x_up-x_low +1);
 363
 364   //Now, calculate the center of mass in x-direction
 365   for(Int_t xbin=x_low; xbin <= x_up; xbin++)
 366     {
 367       value = hist->GetXaxis()->GetBinCenter(xbin);
 368       top += value*m[xbin];
 369       butt += m[xbin];
 370     }
 371   x_peak = top/butt;
 372
 373   //printf("FOund x_peak %f\n",x_peak);
 374
 375   //Find the peak in y direction:
 376   Int_t x_l = hist->GetXaxis()->FindBin(x_peak);
 377   if(hist->GetXaxis()->GetBinCenter(x_l) > x_peak)
 378     x_l--;
 379
 380   Int_t x_u = x_l + 1;
 381
 382   if(hist->GetXaxis()->GetBinCenter(x_l) > x_peak || hist->GetXaxis()->GetBinCenter(x_u) <= x_peak)
 383     printf("\nAliL3HoughMaxFinder::FindPeak : Wrong xrange %f %f %f\n\n",hist->GetXaxis()->GetBinCenter(x_l),x_peak,hist->GetXaxis()->GetBinCenter(x_u));
 384
 385   //printf("\nxlow %f xup %f\n",hist->GetXaxis()->GetBinCenter(x_l),hist->GetXaxis()->GetBinCenter(x_u));
 386
 387   value=top=butt=0;
 388
 389   //printf("ylow %f yup %f\n",hist->GetYaxis()->GetBinCenter(m_low[x_l]),hist->GetYaxis()->GetBinCenter(m_up[x_l]));
 390   //printf("ylow %f yup %f\n",hist->GetYaxis()->GetBinCenter(m_low[x_u]),hist->GetYaxis()->GetBinCenter(m_up[x_u]));
 391
 392   for(Int_t ybin=m_low[x_l]; ybin <= m_up[x_l]; ybin++)
 393     {
 394       value = hist->GetYaxis()->GetBinCenter(ybin);
 395       bin = hist->GetBin(x_l,ybin);
 396       top += value*hist->GetBinContent(bin);
 397       butt += hist->GetBinContent(bin);
 398     }
 399   Double_t y_peak_low = top/butt;
 400
 401   //printf("y_peak_low %f\n",y_peak_low);
 402
 403   value=top=butt=0;
 404   for(Int_t ybin=m_low[x_u]; ybin <= m_up[x_u]; ybin++)
 405     {
 406       value = hist->GetYaxis()->GetBinCenter(ybin);
 407       bin = hist->GetBin(x_u,ybin);
 408       top += value*hist->GetBinContent(bin);
 409       butt += hist->GetBinContent(bin);
 410     }
 411   Double_t y_peak_up = top/butt;
 412
 413   //printf("y_peak_up %f\n",y_peak_up);
 414
 415   Double_t x_value_up = hist->GetXaxis()->GetBinCenter(x_u);
 416   Double_t x_value_low = hist->GetXaxis()->GetBinCenter(x_l);
 417
 418   Double_t y_peak = (y_peak_low*(x_value_up - x_peak) + y_peak_up*(x_peak - x_value_low))/(x_value_up - x_value_low);
 419
 420   //Find the weight:
 421   bin = hist->FindBin(x_peak,y_peak);
 422   Int_t weight = (Int_t)hist->GetBinContent(bin);
 423
 424   if(weight==0) return false;
 425
 426   track->SetTrackParameters(x_peak,y_peak,weight);
 427
 428   delete [] m;
 429   delete [] m_up;
 430   delete [] m_low;
 431   return true;
 432 }
 433
 434
 435 AliL3TrackArray *AliL3HoughMaxFinder::FindPeak(TH2F *hist,Int_t t1,Double_t t2,Int_t t3)
 436 {
 437   //Attempt of a more sophisticated peak finder.
 438
 439   AliL3TrackArray *tracks = new AliL3TrackArray("AliL3HoughTrack",1);
 440
 441   Int_t xmin = hist->GetXaxis()->GetFirst();
 442   Int_t xmax = hist->GetXaxis()->GetLast();
 443   Int_t ymin = hist->GetYaxis()->GetFirst();
 444   Int_t ymax = hist->GetYaxis()->GetLast();
 445   Int_t nbinsx = hist->GetXaxis()->GetNbins()+1;
 446
 447   Int_t *m = new Int_t[nbinsx];
 448   Int_t *m_low = new Int_t[nbinsx];
 449   Int_t *m_up = new Int_t[nbinsx];
 450
 451
 452  recompute:  //this is a goto.
 453
 454   for(Int_t i=0; i<nbinsx; i++)
 455     {
 456       m[i]=0;
 457       m_low[i]=0;
 458       m_up[i]=0;
 459     }
 460
 461   Int_t max_x=0,sum=0,max_xbin=0,bin;
 462
 463   for(Int_t xbin=xmin; xbin<=xmax; xbin++)
 464     {
 465       for(Int_t ybin=ymin; ybin < ymax - t1; ybin++)
 466         {
 467           sum = 0;
 468           for(Int_t y=ybin; y < ybin+t1; y++)
 469             {
 470               //Inside window
 471               bin = hist->GetBin(xbin,y);
 472               sum += (Int_t)hist->GetBinContent(bin);
 473
 474             }
 475           if(sum > m[xbin]) //Max value locally in this xbin
 476             {
 477               m[xbin]=sum;
 478               m_low[xbin]=ybin;
 479               m_up[xbin]=ybin + t1 - 1;
 480             }
 481
 482         }
 483
 484       if(m[xbin] > max_x) //Max value globally in x-direction
 485         {
 486           max_xbin = xbin;
 487           max_x = m[xbin];//sum;
 488         }
 489     }
 490   //printf("max_xbin %d max_x %d m_low %d m_up %d\n",max_xbin,max_x,m_low[max_xbin],m_up[max_xbin]);
 491   //printf("ylow %f yup %f\n",hist->GetYaxis()->GetBinCenter(m_low[max_xbin]),hist->GetYaxis()->GetBinCenter(m_up[max_xbin]));
 492
 493   //Determine a width in the x-direction
 494   Int_t x_low=0,x_up=0;
 495   for(Int_t xbin=max_xbin-1; xbin >= xmin; xbin--)
 496     {
 497       if(m[xbin]*t2 < max_x)
 498         {
 499           x_low = xbin+1;
 500           break;
 501         }
 502     }
 503   for(Int_t xbin = max_xbin+1; xbin <=xmax; xbin++)
 504     {
 505       if(m[xbin]*t2 < max_x)
 506         {
 507           x_up = xbin-1;
 508           break;
 509         }
 510     }
 511
 512   Double_t top=0,butt=0,value,x_peak;
 513   if(x_up - x_low + 1 > t3)
 514     {
 515       t1 -= 1;
 516       printf("\nxrange out if limit x_up %d x_low %d\n\n",x_low,x_up);
 517       if(t1 > 1)
 518         goto recompute;
 519       else
 520         {
 521           x_peak = hist->GetXaxis()->GetBinCenter(max_xbin);
 522           goto moveon;
 523         }
 524     }
 525
 526   //printf("xlow %f xup %f\n",hist->GetXaxis()->GetBinCenter(x_low),hist->GetXaxis()->GetBinCenter(x_up));
 527   //printf("Spread in x %d\n",x_up-x_low +1);
 528
 529   //Now, calculate the center of mass in x-direction
 530   for(Int_t xbin=x_low; xbin <= x_up; xbin++)
 531     {
 532       value = hist->GetXaxis()->GetBinCenter(xbin);
 533       top += value*m[xbin];
 534       butt += m[xbin];
 535     }
 536   x_peak = top/butt;
 537
 538  moveon:
 539
 540   //Find the peak in y direction:
 541   Int_t x_l = hist->GetXaxis()->FindBin(x_peak);
 542   if(hist->GetXaxis()->GetBinCenter(x_l) > x_peak)
 543     x_l--;
 544
 545   Int_t x_u = x_l + 1;
 546
 547   if(hist->GetXaxis()->GetBinCenter(x_l) > x_peak || hist->GetXaxis()->GetBinCenter(x_u) <= x_peak)
 548     printf("\nAliL3HoughMaxFinder::FindPeak : Wrong xrange %f %f %f\n\n",hist->GetXaxis()->GetBinCenter(x_l),x_peak,hist->GetXaxis()->GetBinCenter(x_u));
 549
 550     //printf("\nxlow %f xup %f\n",hist->GetXaxis()->GetBinCenter(x_l),hist->GetXaxis()->GetBinCenter(x_u));
 551
 552   value=top=butt=0;
 553
 554   //printf("ylow %f yup %f\n",hist->GetYaxis()->GetBinCenter(m_low[x_l]),hist->GetYaxis()->GetBinCenter(m_up[x_l]));
 555   //printf("ylow %f yup %f\n",hist->GetYaxis()->GetBinCenter(m_low[x_u]),hist->GetYaxis()->GetBinCenter(m_up[x_u]));
 556
 557   for(Int_t ybin=m_low[x_l]; ybin <= m_up[x_l]; ybin++)
 558     {
 559       value = hist->GetYaxis()->GetBinCenter(ybin);
 560       bin = hist->GetBin(x_l,ybin);
 561       top += value*hist->GetBinContent(bin);
 562       butt += hist->GetBinContent(bin);
 563     }
 564   Double_t y_peak_low = top/butt;
 565
 566   //printf("y_peak_low %f\n",y_peak_low);
 567
 568   value=top=butt=0;
 569   for(Int_t ybin=m_low[x_u]; ybin <= m_up[x_u]; ybin++)
 570     {
 571       value = hist->GetYaxis()->GetBinCenter(ybin);
 572       bin = hist->GetBin(x_u,ybin);
 573       top += value*hist->GetBinContent(bin);
 574       butt += hist->GetBinContent(bin);
 575     }
 576   Double_t y_peak_up = top/butt;
 577
 578   //printf("y_peak_up %f\n",y_peak_up);
 579
 580   Double_t x_value_up = hist->GetXaxis()->GetBinCenter(x_u);
 581   Double_t x_value_low = hist->GetXaxis()->GetBinCenter(x_l);
 582
 583   Double_t y_peak = (y_peak_low*(x_value_up - x_peak) + y_peak_up*(x_peak - x_value_low))/(x_value_up - x_value_low);
 584
 585
 586   //Find the weight:
 587   bin = hist->FindBin(x_peak,y_peak);
 588   Int_t weight = (Int_t)hist->GetBinContent(bin);
 589
 590   AliL3HoughTrack *track = (AliL3HoughTrack*)tracks->NextTrack();
 591   track->SetTrackParameters(x_peak,y_peak,weight);
 592
 593
 594   delete [] m;
 595   delete [] m_low;
 596   delete [] m_up;
 597
 598   return tracks;
 599
 600
 601 }
 602