ROOT logo
/************************************************************************* 
* Copyright(c) 1998-2048, ALICE Experiment at CERN, All rights reserved. * 
*                                                                        * 
* Author: The ALICE Off-line Project.                                    * 
* Contributors are mentioned in the code where appropriate.              * 
*                                                                        * 
* 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.                  * 
**************************************************************************/

// Author: Dariusz Miskowiec <mailto:d.miskowiec@gsi.de> 2007

//=============================================================================
// multidimensional histogram 
// Physically, the data are kept in one single one-dimensional histogram. 
// Projecting on 1, 2, and n-1 dimensions is implemented.
// The histogram can be saved on root file as the one-dimensional data 
// histogram and the axes, thus eternal forward compatibility is ensured. 
//=============================================================================

#include <cmath>
#include <stdlib.h>
#include <TFile.h>
#include <TDirectory.h>
#include <TROOT.h>
#include <TAxis.h>
#include <TH2.h>
#include <TObjArray.h>
#include <TObjString.h>
#include <TString.h>
#include "AliUnicorHN.h"

ClassImp(AliUnicorHN)

//=============================================================================
AliUnicorHN::AliUnicorHN(const char *nam, Int_t ndim, TAxis **ax) 
  : TH1D(nam, nam, Albins(ndim,ax), 0.5, Albins(ndim,ax)+0.5), 
    fNdim(ndim) {

  // constructor
 
  // Above, we just have managed to create a single one-dimensional histogram 
  // with number of bins equal to the product of the numbers of bins in all 
  // dimensions. For easy inspection the histogram range was set to -0.5,n-0.5.

  for (int i=0; i<fNdim; i++) ax[i]->Copy(fAxis[i]); 
  for (int i=0; i<fNdim; i++) fAxis[i].SetName(Form("axis%d",i));
  
  // for speed reasons, number of bins of each axis is stored on fNbins
  // and the running product of the latter on fMbins
  // so fMbins = {...,fNbins[fNdim-2]*fNbins[fNdim-1],fNbins[fNdim-1],1}
  
  for (int i=0; i<fgkMaxNdim; i++) fNbins[i] = fAxis[i].GetNbins();
  for (int i=0; i<fgkMaxNdim; i++) fMbins[i] = 1;
  for (int i=fNdim-1; i>0; i--) fMbins[i-1] = fMbins[i]*fNbins[i];
  printf("   %d-dimensional histogram %s with %d bins created\n",fNdim,nam,GetNbinsX());
}
//=============================================================================
AliUnicorHN* AliUnicorHN::Retrieve(const char *filnam, const char *nam) { 

  // retrieve a multidimensional histogram from file

  TFile *f = TFile::Open(filnam,"read");
  if (!f) return 0;
  if (!f->cd(nam)) {f->Close(); return 0;}
  TH1D *hist = (TH1D*) gDirectory->Get("histo");
  if (!hist) {printf("No histogram histo on file %s in directory %s\n",filnam,nam); return 0;}
  hist->SetDirectory(gROOT);
  TAxis *ax[fgkMaxNdim]={0};
  int n=0;
  while ((ax[n] = (TAxis *) gDirectory->Get(Form("axis%d",n)))) n++; 
  f->Close();
  
  if (hist->GetNbinsX()!=Albins(n,ax)) {
    printf("number of bins of histo %d differs from product of nbins of axes %d\n",
	   hist->GetNbinsX(),Albins(n,ax));
    return 0;
  }

  // derive the name from the path (ccc from aaa/bbb/ccc)
  TString path=nam;
  TObjArray *ar = path.Tokenize("/");
  TObjString *os = (TObjString*) ar->Last();
  const char *lastnam = os->GetString().Data();

  AliUnicorHN *hi = new AliUnicorHN(lastnam,n,ax);
  //  *((TH1D*) hi) = *hist;
  hist->Copy(*((TH1D*)hi));
  hi->SetName(lastnam);
  delete hist;
  delete ar;
  return hi;
}
//=============================================================================
Int_t AliUnicorHN::Albins(Int_t n, TAxis **ax) {

  // Calculate product of nbins of ax[0]...ax[n-1]
  // (= total number of bins of the multidimensional histogram to be made). 

  Int_t nbins = 1;
  //  while (n--) nbins *= ax[n]->GetNbins();
  for (int i=0; i<n; i++) nbins *= ax[i]->GetNbins();
  return nbins;
}
//=============================================================================
Int_t AliUnicorHN::MulToOne(const Int_t * const k) const {

  // Calculate the 1-dim index n from n-dim indices k[fNdim].
  // Valid k[i] should be between 0 and fNbins[i]-1.
  // Valid result will be between 0 and GetNbinsX()-1. 
  // Return -1 if under- or overflow in any dimension.

  Int_t n = 0;
  for (int i=0; i<fNdim; i++) {
    if (k[i]<0) return -1;
    if (k[i]>=fNbins[i]) return -1;
    n += fMbins[i]*k[i];
  }
  return n;
}
//=============================================================================
Int_t AliUnicorHN::MulToOne(Double_t *x) {

  // Calculate the 1-dim index n from n-dim vector x, representing the 
  // abscissa of the n-dim histogram. The result will be between 0 and 
  // GetNbinsX()-1. 

  Int_t k[fgkMaxNdim] = {0};
  for (int i=0; i<fNdim; i++) k[i] = fAxis[i].FindFixBin(x[i])-1;
  return MulToOne(k);
}
//=============================================================================
void AliUnicorHN::OneToMul(Int_t n, Int_t *k) const {

  // Calculate the n-dim indices k[fNdim] from 1-dim index n.
  // Valid n should be between 0 and GetNbinsX()-1. 
  // Valid results will be between 0 and fNbins[i]-1.

  div_t idi; // integer division structure 
  for (int i=0; i<fNdim; i++) {
    idi  = div(n,fMbins[i]);
    k[i] = idi.quot;  // quotient
    n    = idi.rem;   // remainder
  }
}
//=============================================================================
Int_t AliUnicorHN::Fill(Double_t *xx, Double_t w) {

  // Fill the histogram. The array xx holds the abscissa information, w is the 
  // weigth. The 1-dim histogram is filled using the standard Fill method 
  // (direct access to the arrays was tried and was not faster). 

  int nbin = MulToOne(xx);
  if (nbin == -1) return 0;
  return TH1D::Fill(nbin+1,w); 
}
//=============================================================================
Int_t AliUnicorHN::Fill(Double_t x0, Double_t x1, ...) {

  // Fill the histogram. Arguments are passed as doubles rather than array. 

  va_list ap;
  Double_t xx[fgkMaxNdim] = {x0, x1};
  va_start(ap,x1);
  for (int i=2; i<fNdim; i++) xx[i] = va_arg(ap,Double_t);
  Double_t weigth = va_arg(ap,Double_t);
  va_end(ap);
  return Fill(xx,weigth);
}
//=============================================================================
Int_t AliUnicorHN::Save() const {

  // Save the 1-dim histo and the axes in a subdirectory on file. This might 
  // not be the most elegant way but it is very simple and backward and forward 
  // compatible. 

  Int_t nbytes = 0;
  TH1D histo(*this);
  histo.SetName("histo"); // hadd bug fix; actually, does not cost memory, strange

  TDirectory *dest = gDirectory->mkdir(GetName());
  if (dest) {
    dest->cd();
    nbytes += histo.Write("histo");
    for (int i=0; i<fNdim; i++) nbytes += fAxis[i].Write();
    printf("   %s stored in %s\n",GetName(),dest->GetPath());
    dest->cd("..");
  }
  return nbytes;
}
//=============================================================================
AliUnicorHN *AliUnicorHN::ProjectAlong(const char *nam, Int_t dim, Int_t first, Int_t last) {

  // Reduce dimension dim by summing up its bins between first and last. 
  // Use root convention: bin=1 is the first bin, bin=nbins is the last. 
  // last=0 means till the last bin
  // Return the resulting fNdim-1 dimensional histogram. 

  if (dim<0 || dim>fNdim-1) return 0;
  if (last<=0) last = fNbins[dim];

  // create new (reduced) histogram

  TAxis *ax[fgkMaxNdim-1] = {0};
  int n=0;
  for (int i=0; i<fNdim; i++) if (i!=dim) ax[n++] = GetAxis(i);
  const char *name = strlen(nam)? nam : Form("%s_wo%d",GetName(),dim);
  const char *eame = Form("%s_error",name); 
  AliUnicorHN *his = new AliUnicorHN(name,fNdim-1,ax); // result histogram
  AliUnicorHN *eis = new AliUnicorHN(eame,fNdim-1,ax); // temporary storage for errors

  // sum up the content and errors squared

  int k[fgkMaxNdim] = {0}; // old hist multiindex
  int m[fgkMaxNdim] = {0}; // new hist multiindex
  for (int i=0; i<GetNbinsX(); i++) {
    OneToMul(i,k);
    if (k[dim]+1<first) continue;
    if (k[dim]+1>last) continue;
    n = 0;
    for (int j=0; j<fNdim; j++) if (j!=dim) m[n++] = k[j];
    n = his->MulToOne(m);
    his->AddBinContent(n+1,GetBinContent(i+1));
    eis->AddBinContent(n+1,GetBinError(i+1)*GetBinError(i+1));
  }

  // combine content and errors in one histogram

  for (int i=0; i<his->GetNbinsX(); i++) {
    his->SetBinError(i+1,sqrt(eis->GetBinContent(i+1)));
  }

  his->SetLineColor(this->GetLineColor());
  his->SetFillColor(this->GetFillColor());
  his->SetMarkerColor(this->GetMarkerColor());
  his->SetMarkerStyle(this->GetMarkerStyle());

  // some cleanup

  delete eis;
  return his;
}
//=============================================================================
TH1D *AliUnicorHN::ProjectOn(const char *nam, Int_t dim, const Int_t * const first, 
		     const Int_t * const last) const {

  // Project on dimension dim. Use only bins between first[i] and last[i]. 
  // Use root convention: bin=1 is the first bin, bin=nbins is the last. 
  // first[i]=0 means from the first bin
  // last[i]=0 means till the last bin

  if (dim<0 || dim>fNdim-1) return 0;

  // calculate the projection; lowest index 0

  double *yy = new double[fNbins[dim]]; // value
  double *ey = new double[fNbins[dim]]; // error
  for (int i=0; i<fNbins[dim]; i++) yy[i]=0;
  for (int i=0; i<fNbins[dim]; i++) ey[i]=0;
  Int_t *k = new Int_t[fNdim];
  for (int i=0; i<GetNbinsX(); i++) {
    OneToMul(i,k);
    int isgood = 1; // bin within the range?
    for (int j=0; j<fNdim; j++) {
      if (first) if (first[j]) if (k[j]+1<first[j]) {isgood = 0; break;}
      if (last)  if (last[j])  if (k[j]+1>last[j])  {isgood = 0; break;}
    }
    if (isgood) {
      yy[k[dim]]+=GetBinContent(i+1);
      ey[k[dim]]+=GetBinError(i+1)*GetBinError(i+1);
    }
  }

  // make the projection histogram
  // use name nam if specified; otherwise generate one

  TH1D *his;
  const char *name = strlen(nam)? nam : Form("%s_proj%d",GetName(),dim);
  if (fAxis[dim].IsVariableBinSize()) 
    his = new TH1D(name,name,fNbins[dim],fAxis[dim].GetXbins()->GetArray());
  else 
    his = new TH1D(name,name,fNbins[dim],fAxis[dim].GetXmin(),fAxis[dim].GetXmax());
  his->SetXTitle(fAxis[dim].GetTitle());
  // or his->GetXaxis()->ImportAttributes(ax);
  his->Sumw2();
  his->SetLineColor(this->GetLineColor());
  his->SetFillColor(this->GetFillColor());
  his->SetMarkerColor(this->GetMarkerColor());
  his->SetMarkerStyle(this->GetMarkerStyle());
  for (int i=0; i<his->GetNbinsX(); i++) {
    his->SetBinContent(i+1,yy[i]);
    his->SetBinError(i+1,sqrt(ey[i]));
  }

  // some cleanup

  delete [] yy;
  delete [] ey;
  delete [] k;
  //  if (name!=nam) delete [] name;

  return his;
}
//=============================================================================
TH1D *AliUnicorHN::ProjectOn(const char *nam, Int_t dim, const Double_t * const first, 
		     const Double_t * const last) const {

  // Project on dimension dim. Use only bins between first[i] and last[i]. 

  if (dim<0 || dim>fNdim-1) return 0;
  Int_t kfirst[fgkMaxNdim] = {0};
  Int_t klast[fgkMaxNdim] = {0};
  for (int i=0; i<fNdim; i++) {
    kfirst[i] = fAxis[i].FindFixBin(first[i]);
    klast[i] = fAxis[i].FindFixBin(last[i]);
  }
  return ProjectOn(nam,dim,kfirst,klast);
}
//=============================================================================
TH2D *AliUnicorHN::ProjectOn(const char *nam, Int_t dim0, Int_t dim1, const Int_t * 
		     const first, const Int_t * const last) const {

  // Project on dim1 vs dim0. Use only bins between first[i] and last[i]. 
  // Use root convention: bin=1 is the first bin, bin=nbins is the last. 
  // first[i]=0 means from the first bin
  // last[i]=0 means till the last bin

  if (dim0<0 || dim0>fNdim-1) return 0;
  if (dim1<0 || dim1>fNdim-1) return 0;

  // calculate the projection

  double **yy = new double*[fNbins[dim0]]; // value
  double **ey = new double*[fNbins[dim0]]; // error
  for (int i=0; i<fNbins[dim0]; i++) yy[i] = new double[fNbins[dim1]]; 
  for (int i=0; i<fNbins[dim0]; i++) ey[i] = new double[fNbins[dim1]]; 
  for (int i=0; i<fNbins[dim0]; i++) for (int j=0; j<fNbins[dim1]; j++) yy[i][j]=0;
  for (int i=0; i<fNbins[dim0]; i++) for (int j=0; j<fNbins[dim1]; j++) ey[i][j]=0;
  Int_t *k = new Int_t[fNdim];
  for (int i=0; i<GetNbinsX(); i++) {
    OneToMul(i,k);
    int isgood = 1; // bin within the range?
    for (int j=0; j<fNdim; j++) {
      if (first) if (first[j]) if (k[j]+1<first[j]) {isgood = 0; break;}
      if (last)  if (last[j])  if (k[j]+1>last[j])  {isgood = 0; break;}
    }
    if (isgood) {
      yy[k[dim0]][k[dim1]]+=GetBinContent(i+1);
      ey[k[dim0]][k[dim1]]+=GetBinError(i+1)*GetBinError(i+1);
    }
  }

  // make the projection histogram
  // use name nam if specified; otherwise generate one

  TH2D *his=0;
  const char *name = strlen(nam)? nam : Form("%s_proj%dvs%d",GetName(),dim1,dim0);
  if (fAxis[dim0].IsVariableBinSize() && fAxis[dim1].IsVariableBinSize()) 
    his = new TH2D(name,name,
		   fNbins[dim0],fAxis[dim0].GetXbins()->GetArray(),
		   fNbins[dim1],fAxis[dim1].GetXbins()->GetArray());
  else if (!fAxis[dim0].IsVariableBinSize() && fAxis[dim1].IsVariableBinSize()) 
    his = new TH2D(name,name,
		   fNbins[dim0],fAxis[dim0].GetXmin(),fAxis[dim0].GetXmax(),
		   fNbins[dim1],fAxis[dim1].GetXbins()->GetArray());
  else if (fAxis[dim0].IsVariableBinSize() && !fAxis[dim1].IsVariableBinSize()) 
    his = new TH2D(name,name,
		   fNbins[dim0],fAxis[dim0].GetXbins()->GetArray(),
		   fNbins[dim1],fAxis[dim1].GetXmin(),fAxis[dim1].GetXmax());
  else if (!fAxis[dim0].IsVariableBinSize() && !fAxis[dim1].IsVariableBinSize()) 
    his = new TH2D(name,name,
		   fNbins[dim0],fAxis[dim0].GetXmin(),fAxis[dim0].GetXmax(),
		   fNbins[dim1],fAxis[dim1].GetXmin(),fAxis[dim1].GetXmax());
  his->SetXTitle(fAxis[dim0].GetTitle());
  his->SetYTitle(fAxis[dim1].GetTitle());
  // or his->GetXaxis()->ImportAttributes(ax0);
  // or his->GetYaxis()->ImportAttributes(ax1);
  his->Sumw2();
  his->SetLineColor(this->GetLineColor());
  his->SetFillColor(this->GetFillColor());
  his->SetMarkerColor(this->GetMarkerColor());
  his->SetMarkerStyle(this->GetMarkerStyle());
  for (int i=0; i<his->GetNbinsX(); i++) 
  for (int j=0; j<his->GetNbinsY(); j++) {
    his->SetBinContent(i+1,j+1,yy[i][j]);
    his->SetBinError(i+1,j+1,sqrt(ey[i][j]));
  }

  // some cleanup

  for (int i=0; i<fNbins[dim0]; i++) delete [] yy[i];
  for (int i=0; i<fNbins[dim0]; i++) delete [] ey[i];
  delete [] yy;
  delete [] ey;
  delete [] k;
  //  if (name!=nam) delete [] name;

  return his;
}
//=============================================================================
TH2D *AliUnicorHN::ProjectOn(const char *nam, Int_t dim0, Int_t dim1, const Double_t * 
		     const first, const Double_t * const last) const {

  // Project on dim1 vs dim0. Use only bins between first[i] and last[i]. 

  if (dim0<0 || dim0>fNdim-1) return 0;
  if (dim1<0 || dim1>fNdim-1) return 0;
  Int_t kfirst[fgkMaxNdim] = {0};
  Int_t klast[fgkMaxNdim] = {0};
  for (int i=0; i<fNdim; i++) {
    kfirst[i] = fAxis[i].FindFixBin(first[i]);
    klast[i] = fAxis[i].FindFixBin(last[i]);
  }
  return ProjectOn(nam,dim0,dim1,kfirst,klast);
}
//=============================================================================
 AliUnicorHN.cxx:1
 AliUnicorHN.cxx:2
 AliUnicorHN.cxx:3
 AliUnicorHN.cxx:4
 AliUnicorHN.cxx:5
 AliUnicorHN.cxx:6
 AliUnicorHN.cxx:7
 AliUnicorHN.cxx:8
 AliUnicorHN.cxx:9
 AliUnicorHN.cxx:10
 AliUnicorHN.cxx:11
 AliUnicorHN.cxx:12
 AliUnicorHN.cxx:13
 AliUnicorHN.cxx:14
 AliUnicorHN.cxx:15
 AliUnicorHN.cxx:16
 AliUnicorHN.cxx:17
 AliUnicorHN.cxx:18
 AliUnicorHN.cxx:19
 AliUnicorHN.cxx:20
 AliUnicorHN.cxx:21
 AliUnicorHN.cxx:22
 AliUnicorHN.cxx:23
 AliUnicorHN.cxx:24
 AliUnicorHN.cxx:25
 AliUnicorHN.cxx:26
 AliUnicorHN.cxx:27
 AliUnicorHN.cxx:28
 AliUnicorHN.cxx:29
 AliUnicorHN.cxx:30
 AliUnicorHN.cxx:31
 AliUnicorHN.cxx:32
 AliUnicorHN.cxx:33
 AliUnicorHN.cxx:34
 AliUnicorHN.cxx:35
 AliUnicorHN.cxx:36
 AliUnicorHN.cxx:37
 AliUnicorHN.cxx:38
 AliUnicorHN.cxx:39
 AliUnicorHN.cxx:40
 AliUnicorHN.cxx:41
 AliUnicorHN.cxx:42
 AliUnicorHN.cxx:43
 AliUnicorHN.cxx:44
 AliUnicorHN.cxx:45
 AliUnicorHN.cxx:46
 AliUnicorHN.cxx:47
 AliUnicorHN.cxx:48
 AliUnicorHN.cxx:49
 AliUnicorHN.cxx:50
 AliUnicorHN.cxx:51
 AliUnicorHN.cxx:52
 AliUnicorHN.cxx:53
 AliUnicorHN.cxx:54
 AliUnicorHN.cxx:55
 AliUnicorHN.cxx:56
 AliUnicorHN.cxx:57
 AliUnicorHN.cxx:58
 AliUnicorHN.cxx:59
 AliUnicorHN.cxx:60
 AliUnicorHN.cxx:61
 AliUnicorHN.cxx:62
 AliUnicorHN.cxx:63
 AliUnicorHN.cxx:64
 AliUnicorHN.cxx:65
 AliUnicorHN.cxx:66
 AliUnicorHN.cxx:67
 AliUnicorHN.cxx:68
 AliUnicorHN.cxx:69
 AliUnicorHN.cxx:70
 AliUnicorHN.cxx:71
 AliUnicorHN.cxx:72
 AliUnicorHN.cxx:73
 AliUnicorHN.cxx:74
 AliUnicorHN.cxx:75
 AliUnicorHN.cxx:76
 AliUnicorHN.cxx:77
 AliUnicorHN.cxx:78
 AliUnicorHN.cxx:79
 AliUnicorHN.cxx:80
 AliUnicorHN.cxx:81
 AliUnicorHN.cxx:82
 AliUnicorHN.cxx:83
 AliUnicorHN.cxx:84
 AliUnicorHN.cxx:85
 AliUnicorHN.cxx:86
 AliUnicorHN.cxx:87
 AliUnicorHN.cxx:88
 AliUnicorHN.cxx:89
 AliUnicorHN.cxx:90
 AliUnicorHN.cxx:91
 AliUnicorHN.cxx:92
 AliUnicorHN.cxx:93
 AliUnicorHN.cxx:94
 AliUnicorHN.cxx:95
 AliUnicorHN.cxx:96
 AliUnicorHN.cxx:97
 AliUnicorHN.cxx:98
 AliUnicorHN.cxx:99
 AliUnicorHN.cxx:100
 AliUnicorHN.cxx:101
 AliUnicorHN.cxx:102
 AliUnicorHN.cxx:103
 AliUnicorHN.cxx:104
 AliUnicorHN.cxx:105
 AliUnicorHN.cxx:106
 AliUnicorHN.cxx:107
 AliUnicorHN.cxx:108
 AliUnicorHN.cxx:109
 AliUnicorHN.cxx:110
 AliUnicorHN.cxx:111
 AliUnicorHN.cxx:112
 AliUnicorHN.cxx:113
 AliUnicorHN.cxx:114
 AliUnicorHN.cxx:115
 AliUnicorHN.cxx:116
 AliUnicorHN.cxx:117
 AliUnicorHN.cxx:118
 AliUnicorHN.cxx:119
 AliUnicorHN.cxx:120
 AliUnicorHN.cxx:121
 AliUnicorHN.cxx:122
 AliUnicorHN.cxx:123
 AliUnicorHN.cxx:124
 AliUnicorHN.cxx:125
 AliUnicorHN.cxx:126
 AliUnicorHN.cxx:127
 AliUnicorHN.cxx:128
 AliUnicorHN.cxx:129
 AliUnicorHN.cxx:130
 AliUnicorHN.cxx:131
 AliUnicorHN.cxx:132
 AliUnicorHN.cxx:133
 AliUnicorHN.cxx:134
 AliUnicorHN.cxx:135
 AliUnicorHN.cxx:136
 AliUnicorHN.cxx:137
 AliUnicorHN.cxx:138
 AliUnicorHN.cxx:139
 AliUnicorHN.cxx:140
 AliUnicorHN.cxx:141
 AliUnicorHN.cxx:142
 AliUnicorHN.cxx:143
 AliUnicorHN.cxx:144
 AliUnicorHN.cxx:145
 AliUnicorHN.cxx:146
 AliUnicorHN.cxx:147
 AliUnicorHN.cxx:148
 AliUnicorHN.cxx:149
 AliUnicorHN.cxx:150
 AliUnicorHN.cxx:151
 AliUnicorHN.cxx:152
 AliUnicorHN.cxx:153
 AliUnicorHN.cxx:154
 AliUnicorHN.cxx:155
 AliUnicorHN.cxx:156
 AliUnicorHN.cxx:157
 AliUnicorHN.cxx:158
 AliUnicorHN.cxx:159
 AliUnicorHN.cxx:160
 AliUnicorHN.cxx:161
 AliUnicorHN.cxx:162
 AliUnicorHN.cxx:163
 AliUnicorHN.cxx:164
 AliUnicorHN.cxx:165
 AliUnicorHN.cxx:166
 AliUnicorHN.cxx:167
 AliUnicorHN.cxx:168
 AliUnicorHN.cxx:169
 AliUnicorHN.cxx:170
 AliUnicorHN.cxx:171
 AliUnicorHN.cxx:172
 AliUnicorHN.cxx:173
 AliUnicorHN.cxx:174
 AliUnicorHN.cxx:175
 AliUnicorHN.cxx:176
 AliUnicorHN.cxx:177
 AliUnicorHN.cxx:178
 AliUnicorHN.cxx:179
 AliUnicorHN.cxx:180
 AliUnicorHN.cxx:181
 AliUnicorHN.cxx:182
 AliUnicorHN.cxx:183
 AliUnicorHN.cxx:184
 AliUnicorHN.cxx:185
 AliUnicorHN.cxx:186
 AliUnicorHN.cxx:187
 AliUnicorHN.cxx:188
 AliUnicorHN.cxx:189
 AliUnicorHN.cxx:190
 AliUnicorHN.cxx:191
 AliUnicorHN.cxx:192
 AliUnicorHN.cxx:193
 AliUnicorHN.cxx:194
 AliUnicorHN.cxx:195
 AliUnicorHN.cxx:196
 AliUnicorHN.cxx:197
 AliUnicorHN.cxx:198
 AliUnicorHN.cxx:199
 AliUnicorHN.cxx:200
 AliUnicorHN.cxx:201
 AliUnicorHN.cxx:202
 AliUnicorHN.cxx:203
 AliUnicorHN.cxx:204
 AliUnicorHN.cxx:205
 AliUnicorHN.cxx:206
 AliUnicorHN.cxx:207
 AliUnicorHN.cxx:208
 AliUnicorHN.cxx:209
 AliUnicorHN.cxx:210
 AliUnicorHN.cxx:211
 AliUnicorHN.cxx:212
 AliUnicorHN.cxx:213
 AliUnicorHN.cxx:214
 AliUnicorHN.cxx:215
 AliUnicorHN.cxx:216
 AliUnicorHN.cxx:217
 AliUnicorHN.cxx:218
 AliUnicorHN.cxx:219
 AliUnicorHN.cxx:220
 AliUnicorHN.cxx:221
 AliUnicorHN.cxx:222
 AliUnicorHN.cxx:223
 AliUnicorHN.cxx:224
 AliUnicorHN.cxx:225
 AliUnicorHN.cxx:226
 AliUnicorHN.cxx:227
 AliUnicorHN.cxx:228
 AliUnicorHN.cxx:229
 AliUnicorHN.cxx:230
 AliUnicorHN.cxx:231
 AliUnicorHN.cxx:232
 AliUnicorHN.cxx:233
 AliUnicorHN.cxx:234
 AliUnicorHN.cxx:235
 AliUnicorHN.cxx:236
 AliUnicorHN.cxx:237
 AliUnicorHN.cxx:238
 AliUnicorHN.cxx:239
 AliUnicorHN.cxx:240
 AliUnicorHN.cxx:241
 AliUnicorHN.cxx:242
 AliUnicorHN.cxx:243
 AliUnicorHN.cxx:244
 AliUnicorHN.cxx:245
 AliUnicorHN.cxx:246
 AliUnicorHN.cxx:247
 AliUnicorHN.cxx:248
 AliUnicorHN.cxx:249
 AliUnicorHN.cxx:250
 AliUnicorHN.cxx:251
 AliUnicorHN.cxx:252
 AliUnicorHN.cxx:253
 AliUnicorHN.cxx:254
 AliUnicorHN.cxx:255
 AliUnicorHN.cxx:256
 AliUnicorHN.cxx:257
 AliUnicorHN.cxx:258
 AliUnicorHN.cxx:259
 AliUnicorHN.cxx:260
 AliUnicorHN.cxx:261
 AliUnicorHN.cxx:262
 AliUnicorHN.cxx:263
 AliUnicorHN.cxx:264
 AliUnicorHN.cxx:265
 AliUnicorHN.cxx:266
 AliUnicorHN.cxx:267
 AliUnicorHN.cxx:268
 AliUnicorHN.cxx:269
 AliUnicorHN.cxx:270
 AliUnicorHN.cxx:271
 AliUnicorHN.cxx:272
 AliUnicorHN.cxx:273
 AliUnicorHN.cxx:274
 AliUnicorHN.cxx:275
 AliUnicorHN.cxx:276
 AliUnicorHN.cxx:277
 AliUnicorHN.cxx:278
 AliUnicorHN.cxx:279
 AliUnicorHN.cxx:280
 AliUnicorHN.cxx:281
 AliUnicorHN.cxx:282
 AliUnicorHN.cxx:283
 AliUnicorHN.cxx:284
 AliUnicorHN.cxx:285
 AliUnicorHN.cxx:286
 AliUnicorHN.cxx:287
 AliUnicorHN.cxx:288
 AliUnicorHN.cxx:289
 AliUnicorHN.cxx:290
 AliUnicorHN.cxx:291
 AliUnicorHN.cxx:292
 AliUnicorHN.cxx:293
 AliUnicorHN.cxx:294
 AliUnicorHN.cxx:295
 AliUnicorHN.cxx:296
 AliUnicorHN.cxx:297
 AliUnicorHN.cxx:298
 AliUnicorHN.cxx:299
 AliUnicorHN.cxx:300
 AliUnicorHN.cxx:301
 AliUnicorHN.cxx:302
 AliUnicorHN.cxx:303
 AliUnicorHN.cxx:304
 AliUnicorHN.cxx:305
 AliUnicorHN.cxx:306
 AliUnicorHN.cxx:307
 AliUnicorHN.cxx:308
 AliUnicorHN.cxx:309
 AliUnicorHN.cxx:310
 AliUnicorHN.cxx:311
 AliUnicorHN.cxx:312
 AliUnicorHN.cxx:313
 AliUnicorHN.cxx:314
 AliUnicorHN.cxx:315
 AliUnicorHN.cxx:316
 AliUnicorHN.cxx:317
 AliUnicorHN.cxx:318
 AliUnicorHN.cxx:319
 AliUnicorHN.cxx:320
 AliUnicorHN.cxx:321
 AliUnicorHN.cxx:322
 AliUnicorHN.cxx:323
 AliUnicorHN.cxx:324
 AliUnicorHN.cxx:325
 AliUnicorHN.cxx:326
 AliUnicorHN.cxx:327
 AliUnicorHN.cxx:328
 AliUnicorHN.cxx:329
 AliUnicorHN.cxx:330
 AliUnicorHN.cxx:331
 AliUnicorHN.cxx:332
 AliUnicorHN.cxx:333
 AliUnicorHN.cxx:334
 AliUnicorHN.cxx:335
 AliUnicorHN.cxx:336
 AliUnicorHN.cxx:337
 AliUnicorHN.cxx:338
 AliUnicorHN.cxx:339
 AliUnicorHN.cxx:340
 AliUnicorHN.cxx:341
 AliUnicorHN.cxx:342
 AliUnicorHN.cxx:343
 AliUnicorHN.cxx:344
 AliUnicorHN.cxx:345
 AliUnicorHN.cxx:346
 AliUnicorHN.cxx:347
 AliUnicorHN.cxx:348
 AliUnicorHN.cxx:349
 AliUnicorHN.cxx:350
 AliUnicorHN.cxx:351
 AliUnicorHN.cxx:352
 AliUnicorHN.cxx:353
 AliUnicorHN.cxx:354
 AliUnicorHN.cxx:355
 AliUnicorHN.cxx:356
 AliUnicorHN.cxx:357
 AliUnicorHN.cxx:358
 AliUnicorHN.cxx:359
 AliUnicorHN.cxx:360
 AliUnicorHN.cxx:361
 AliUnicorHN.cxx:362
 AliUnicorHN.cxx:363
 AliUnicorHN.cxx:364
 AliUnicorHN.cxx:365
 AliUnicorHN.cxx:366
 AliUnicorHN.cxx:367
 AliUnicorHN.cxx:368
 AliUnicorHN.cxx:369
 AliUnicorHN.cxx:370
 AliUnicorHN.cxx:371
 AliUnicorHN.cxx:372
 AliUnicorHN.cxx:373
 AliUnicorHN.cxx:374
 AliUnicorHN.cxx:375
 AliUnicorHN.cxx:376
 AliUnicorHN.cxx:377
 AliUnicorHN.cxx:378
 AliUnicorHN.cxx:379
 AliUnicorHN.cxx:380
 AliUnicorHN.cxx:381
 AliUnicorHN.cxx:382
 AliUnicorHN.cxx:383
 AliUnicorHN.cxx:384
 AliUnicorHN.cxx:385
 AliUnicorHN.cxx:386
 AliUnicorHN.cxx:387
 AliUnicorHN.cxx:388
 AliUnicorHN.cxx:389
 AliUnicorHN.cxx:390
 AliUnicorHN.cxx:391
 AliUnicorHN.cxx:392
 AliUnicorHN.cxx:393
 AliUnicorHN.cxx:394
 AliUnicorHN.cxx:395
 AliUnicorHN.cxx:396
 AliUnicorHN.cxx:397
 AliUnicorHN.cxx:398
 AliUnicorHN.cxx:399
 AliUnicorHN.cxx:400
 AliUnicorHN.cxx:401
 AliUnicorHN.cxx:402
 AliUnicorHN.cxx:403
 AliUnicorHN.cxx:404
 AliUnicorHN.cxx:405
 AliUnicorHN.cxx:406
 AliUnicorHN.cxx:407
 AliUnicorHN.cxx:408
 AliUnicorHN.cxx:409
 AliUnicorHN.cxx:410
 AliUnicorHN.cxx:411
 AliUnicorHN.cxx:412
 AliUnicorHN.cxx:413
 AliUnicorHN.cxx:414
 AliUnicorHN.cxx:415
 AliUnicorHN.cxx:416
 AliUnicorHN.cxx:417
 AliUnicorHN.cxx:418
 AliUnicorHN.cxx:419
 AliUnicorHN.cxx:420
 AliUnicorHN.cxx:421