ROOT logo
////////////////////////////////////////////////////////
//
// Bucket representation for TKDInterpolator(Base)
//
// The class store data and provides the interface to draw and print.
// The bucket - generalized histogram bin in N dimensions is represented by unprocessed data like
//   - experimental PDF value and statistical error 
//   - COG position (n-tuplu)
//   - boundaries
// and interpolated data like
//   - parameters of the local parabolic fit
//   - their covariance matrix
//  
// For drawing 2D projections the helper class TKDNodeInfo::TKDNodeDraw is used.
//
// Author Alexandru Bercuci <A.Bercuci@gsi.de>
//
////////////////////////////////////////////////////////

#include "TKDNodeInfo.h"

#include "TVectorD.h"
#include "TMatrixD.h"
#include "TRandom.h"
#include "TMath.h"

ClassImp(TKDNodeInfo)
ClassImp(TKDNodeInfo::TKDNodeDraw)


//_________________________________________________________________
TKDNodeInfo::TKDNodeInfo(Int_t dim):
  TObject()
  ,fNDim(3*dim)
  ,fData(NULL)
  ,fNpar(0)
  ,fNcov(0)
  ,fPar(NULL)
  ,fCov(NULL)
{
  // Default constructor
  fVal[0] = 0.; fVal[1] = 0.;
  Build(dim);
}

//_________________________________________________________________
TKDNodeInfo::TKDNodeInfo(const TKDNodeInfo &ref):
  TObject((TObject&) ref)
  ,fNDim(ref.fNDim)
  ,fData(NULL)
  ,fNpar(0)
  ,fNcov(0)
  ,fPar(NULL)
  ,fCov(NULL)
{
  // Copy constructor
  Build(fNDim/3);

  fData = new Float_t[fNDim];
  memcpy(fData, ref.fData, fNDim*sizeof(Float_t));
  fVal[0] = ref.fVal[0];
  fVal[1] = ref.fVal[1];
  if(ref.fNpar&&ref.fPar){ 
    fNpar = ref.fNpar;
    fPar=new Double_t[fNpar];
    memcpy(fPar, ref.fPar, fNpar*sizeof(Double_t));
  }
  if(ref.fNcov && ref.fCov){ 
    fNcov = ref.fNcov;
    fCov=new Double_t[fNcov];
    memcpy(fCov, ref.fCov, fNcov*sizeof(Double_t));
  }
}


//_________________________________________________________________
TKDNodeInfo::~TKDNodeInfo()
{
  // Destructor
  if(fData) delete [] fData;
  if(fPar) delete [] fPar;
  if(fCov) delete [] fCov;
}

//_________________________________________________________________
TKDNodeInfo& TKDNodeInfo::operator=(const TKDNodeInfo & ref)
{
//	Info("operator==()", "...");
  
  if(this == &ref) return *this;
  Int_t ndim = fNDim/3;
  if(fNDim != ref.fNDim){
    fNDim = ref.fNDim;
    Build(ndim);
  }
  memcpy(fData, ref.fData, fNDim*sizeof(Float_t));
  fVal[0] = ref.fVal[0];
  fVal[1] = ref.fVal[1];
  if(ref.fNpar&&ref.fPar){ 
    fNpar = ref.fNpar;
    fPar=new Double_t[fNpar];
    memcpy(fPar, ref.fPar, fNpar*sizeof(Double_t));
  }
  if(ref.fNcov && ref.fCov){ 
    fNcov = ref.fNcov;
    fCov=new Double_t[fNcov];
    memcpy(fCov, ref.fCov, fNcov*sizeof(Double_t));
  }
  return *this;
}

//_________________________________________________________________
void TKDNodeInfo::Build(Int_t dim)
{
// Allocate/Reallocate space for this node.

//	Info("Build()", "...");

  if(!dim) return;
  fNDim = 3*dim;
  if(fData) delete [] fData;
  fData = new Float_t[fNDim];
  return;
}

//_________________________________________________________________
void TKDNodeInfo::Bootstrap()
{
  if(!fNpar || !fPar) return;

  Int_t ndim = Int_t(.5*(TMath::Sqrt(1.+8.*fNpar)-1.))-1;
  fNDim = 3*ndim;
}

//_________________________________________________________________
void TKDNodeInfo::SetNode(Int_t ndim, Float_t *data, Float_t *pdf)
{
  Build(ndim);
  memcpy(fData, data, fNDim*sizeof(Float_t));
  fVal[0]=pdf[0]; fVal[1]=pdf[1];
}


//_________________________________________________________________
void TKDNodeInfo::Print(const Option_t *opt) const
{
  // Print the content of the node
  Int_t dim = Int_t(fNDim/3.);
  printf("x[");
  for(int idim=0; idim<dim; idim++) printf("%f ", fData?fData[idim]:0.);
  printf("] f = [%f +- %f]\n", fVal[0], fVal[1]);

  if(fData){
    Float_t *bounds = &fData[dim];
    printf("range[");
    for(int idim=0; idim<dim; idim++) printf("(%f %f) ", bounds[2*idim], bounds[2*idim+1]);
    printf("]\n");
  }
  if(strcmp(opt, "a")!=0) return;

  if(fNpar){ 
    printf("Fit parameters : \n");
    for(int ip=0; ip<fNpar; ip++) printf("p%d[%f] ", ip, fPar[ip]);
    printf("\n");
  }
  if(!fNcov) return;
  for(int ip(0), n(0); ip<fNpar; ip++){
    for(int jp(ip); jp<fNpar; jp++) printf("c(%d %d)[%f] ", ip, jp, fCov[n++]);
    printf("\n");
  }
}

//_________________________________________________________________
void TKDNodeInfo::Store(TVectorD const *par, TMatrixD const *cov)
{
// Store the parameters and the covariance in the node

  if(!fPar){SetNpar(); fPar = new Double_t[fNpar];}
  for(int ip=0; ip<fNpar; ip++) fPar[ip] = (*par)[ip];

  if(!cov) return;
  if(!fCov){SetNcov(); fCov = new Double_t[fNcov];}
  for(int ip(0), np(0); ip<fNpar; ip++)
    for(int jp=ip; jp<fNpar; jp++) fCov[np++] = (*cov)(ip,jp);
}

//_________________________________________________________________
Bool_t TKDNodeInfo::CookPDF(const Double_t *point, Double_t &result, Double_t &error) const
{
// Recalculate the PDF for one node from the results of interpolation (parameters and covariance matrix)

  Int_t ndim = Int_t(fNDim/3.);
  if(ndim>10) return kFALSE; // support only up to 10 dimensions
  //printf("ndim[%d] npar[%d] ncov[%d]\n", ndim, fNpar, fNcov);

  Double_t fdfdp[66]; memset(fdfdp, 0, ndim*sizeof(Double_t));
  Int_t ipar = 0;
  fdfdp[ipar++] = 1.;
  for(int idim=0; idim<ndim; idim++){
    fdfdp[ipar++] = point[idim];
    for(int jdim=idim; jdim<ndim; jdim++) fdfdp[ipar++] = point[idim]*point[jdim];
  }

  // calculate estimation
  result =0.; error = 0.;
  for(int i=0; i<fNpar; i++) result += fdfdp[i]*fPar[i];
  if(!fNcov) return kTRUE;

  for(int i(0), n(0); i<fNpar; i++){
    error += fdfdp[i]*fdfdp[i]*fCov[n++];
    for(int j(i+1); j<fNpar; j++) error += 2.*fdfdp[i]*fdfdp[j]*fCov[n++];
  }	
  error = TMath::Sqrt(error);
  
  //printf("TKDNodeInfo::CookPDF() : %6.3f +- %6.3f\n", result, error);

  return kTRUE;
}



//_________________________________________________________________
TKDNodeInfo::TKDNodeDraw::TKDNodeDraw() 
  :TBox()
  ,fCOG()
  ,fNode(NULL)
{
  SetFillStyle(3002);
  SetFillColor(50+Int_t(gRandom->Uniform()*50.));

  fCOG.SetMarkerStyle(3);
  fCOG.SetMarkerSize(.7);
  fCOG.SetMarkerColor(2);
}


//_________________________________________________________________
void TKDNodeInfo::TKDNodeDraw::Draw(Option_t* option)
{
  TBox::Draw(option);
  fCOG.Draw("p");
}

//_________________________________________________________________
void TKDNodeInfo::TKDNodeDraw::SetNode(TKDNodeInfo *node, UChar_t size, UChar_t ax1, UChar_t ax2)
{
  fNode=node;
  const Float_t kBorder = 0.;//1.E-4;
  Float_t *bounds = &(node->Data()[size]);
  fX1=bounds[2*ax1]+kBorder;
  fX2=bounds[2*ax1+1]-kBorder;
  fY1=bounds[2*ax2]+kBorder;
  fY2=bounds[2*ax2+1]-kBorder;
  
  Float_t x(node->Data()[ax1]), y(node->Data()[ax2]);
  fCOG.SetX(x); fCOG.SetY(y);
}


//_________________________________________________________________
void TKDNodeInfo::TKDNodeDraw::Print(const Option_t* option) const
{
  if(!fNode) return;
  fNode->Print(option);
}
 TKDNodeInfo.cxx:1
 TKDNodeInfo.cxx:2
 TKDNodeInfo.cxx:3
 TKDNodeInfo.cxx:4
 TKDNodeInfo.cxx:5
 TKDNodeInfo.cxx:6
 TKDNodeInfo.cxx:7
 TKDNodeInfo.cxx:8
 TKDNodeInfo.cxx:9
 TKDNodeInfo.cxx:10
 TKDNodeInfo.cxx:11
 TKDNodeInfo.cxx:12
 TKDNodeInfo.cxx:13
 TKDNodeInfo.cxx:14
 TKDNodeInfo.cxx:15
 TKDNodeInfo.cxx:16
 TKDNodeInfo.cxx:17
 TKDNodeInfo.cxx:18
 TKDNodeInfo.cxx:19
 TKDNodeInfo.cxx:20
 TKDNodeInfo.cxx:21
 TKDNodeInfo.cxx:22
 TKDNodeInfo.cxx:23
 TKDNodeInfo.cxx:24
 TKDNodeInfo.cxx:25
 TKDNodeInfo.cxx:26
 TKDNodeInfo.cxx:27
 TKDNodeInfo.cxx:28
 TKDNodeInfo.cxx:29
 TKDNodeInfo.cxx:30
 TKDNodeInfo.cxx:31
 TKDNodeInfo.cxx:32
 TKDNodeInfo.cxx:33
 TKDNodeInfo.cxx:34
 TKDNodeInfo.cxx:35
 TKDNodeInfo.cxx:36
 TKDNodeInfo.cxx:37
 TKDNodeInfo.cxx:38
 TKDNodeInfo.cxx:39
 TKDNodeInfo.cxx:40
 TKDNodeInfo.cxx:41
 TKDNodeInfo.cxx:42
 TKDNodeInfo.cxx:43
 TKDNodeInfo.cxx:44
 TKDNodeInfo.cxx:45
 TKDNodeInfo.cxx:46
 TKDNodeInfo.cxx:47
 TKDNodeInfo.cxx:48
 TKDNodeInfo.cxx:49
 TKDNodeInfo.cxx:50
 TKDNodeInfo.cxx:51
 TKDNodeInfo.cxx:52
 TKDNodeInfo.cxx:53
 TKDNodeInfo.cxx:54
 TKDNodeInfo.cxx:55
 TKDNodeInfo.cxx:56
 TKDNodeInfo.cxx:57
 TKDNodeInfo.cxx:58
 TKDNodeInfo.cxx:59
 TKDNodeInfo.cxx:60
 TKDNodeInfo.cxx:61
 TKDNodeInfo.cxx:62
 TKDNodeInfo.cxx:63
 TKDNodeInfo.cxx:64
 TKDNodeInfo.cxx:65
 TKDNodeInfo.cxx:66
 TKDNodeInfo.cxx:67
 TKDNodeInfo.cxx:68
 TKDNodeInfo.cxx:69
 TKDNodeInfo.cxx:70
 TKDNodeInfo.cxx:71
 TKDNodeInfo.cxx:72
 TKDNodeInfo.cxx:73
 TKDNodeInfo.cxx:74
 TKDNodeInfo.cxx:75
 TKDNodeInfo.cxx:76
 TKDNodeInfo.cxx:77
 TKDNodeInfo.cxx:78
 TKDNodeInfo.cxx:79
 TKDNodeInfo.cxx:80
 TKDNodeInfo.cxx:81
 TKDNodeInfo.cxx:82
 TKDNodeInfo.cxx:83
 TKDNodeInfo.cxx:84
 TKDNodeInfo.cxx:85
 TKDNodeInfo.cxx:86
 TKDNodeInfo.cxx:87
 TKDNodeInfo.cxx:88
 TKDNodeInfo.cxx:89
 TKDNodeInfo.cxx:90
 TKDNodeInfo.cxx:91
 TKDNodeInfo.cxx:92
 TKDNodeInfo.cxx:93
 TKDNodeInfo.cxx:94
 TKDNodeInfo.cxx:95
 TKDNodeInfo.cxx:96
 TKDNodeInfo.cxx:97
 TKDNodeInfo.cxx:98
 TKDNodeInfo.cxx:99
 TKDNodeInfo.cxx:100
 TKDNodeInfo.cxx:101
 TKDNodeInfo.cxx:102
 TKDNodeInfo.cxx:103
 TKDNodeInfo.cxx:104
 TKDNodeInfo.cxx:105
 TKDNodeInfo.cxx:106
 TKDNodeInfo.cxx:107
 TKDNodeInfo.cxx:108
 TKDNodeInfo.cxx:109
 TKDNodeInfo.cxx:110
 TKDNodeInfo.cxx:111
 TKDNodeInfo.cxx:112
 TKDNodeInfo.cxx:113
 TKDNodeInfo.cxx:114
 TKDNodeInfo.cxx:115
 TKDNodeInfo.cxx:116
 TKDNodeInfo.cxx:117
 TKDNodeInfo.cxx:118
 TKDNodeInfo.cxx:119
 TKDNodeInfo.cxx:120
 TKDNodeInfo.cxx:121
 TKDNodeInfo.cxx:122
 TKDNodeInfo.cxx:123
 TKDNodeInfo.cxx:124
 TKDNodeInfo.cxx:125
 TKDNodeInfo.cxx:126
 TKDNodeInfo.cxx:127
 TKDNodeInfo.cxx:128
 TKDNodeInfo.cxx:129
 TKDNodeInfo.cxx:130
 TKDNodeInfo.cxx:131
 TKDNodeInfo.cxx:132
 TKDNodeInfo.cxx:133
 TKDNodeInfo.cxx:134
 TKDNodeInfo.cxx:135
 TKDNodeInfo.cxx:136
 TKDNodeInfo.cxx:137
 TKDNodeInfo.cxx:138
 TKDNodeInfo.cxx:139
 TKDNodeInfo.cxx:140
 TKDNodeInfo.cxx:141
 TKDNodeInfo.cxx:142
 TKDNodeInfo.cxx:143
 TKDNodeInfo.cxx:144
 TKDNodeInfo.cxx:145
 TKDNodeInfo.cxx:146
 TKDNodeInfo.cxx:147
 TKDNodeInfo.cxx:148
 TKDNodeInfo.cxx:149
 TKDNodeInfo.cxx:150
 TKDNodeInfo.cxx:151
 TKDNodeInfo.cxx:152
 TKDNodeInfo.cxx:153
 TKDNodeInfo.cxx:154
 TKDNodeInfo.cxx:155
 TKDNodeInfo.cxx:156
 TKDNodeInfo.cxx:157
 TKDNodeInfo.cxx:158
 TKDNodeInfo.cxx:159
 TKDNodeInfo.cxx:160
 TKDNodeInfo.cxx:161
 TKDNodeInfo.cxx:162
 TKDNodeInfo.cxx:163
 TKDNodeInfo.cxx:164
 TKDNodeInfo.cxx:165
 TKDNodeInfo.cxx:166
 TKDNodeInfo.cxx:167
 TKDNodeInfo.cxx:168
 TKDNodeInfo.cxx:169
 TKDNodeInfo.cxx:170
 TKDNodeInfo.cxx:171
 TKDNodeInfo.cxx:172
 TKDNodeInfo.cxx:173
 TKDNodeInfo.cxx:174
 TKDNodeInfo.cxx:175
 TKDNodeInfo.cxx:176
 TKDNodeInfo.cxx:177
 TKDNodeInfo.cxx:178
 TKDNodeInfo.cxx:179
 TKDNodeInfo.cxx:180
 TKDNodeInfo.cxx:181
 TKDNodeInfo.cxx:182
 TKDNodeInfo.cxx:183
 TKDNodeInfo.cxx:184
 TKDNodeInfo.cxx:185
 TKDNodeInfo.cxx:186
 TKDNodeInfo.cxx:187
 TKDNodeInfo.cxx:188
 TKDNodeInfo.cxx:189
 TKDNodeInfo.cxx:190
 TKDNodeInfo.cxx:191
 TKDNodeInfo.cxx:192
 TKDNodeInfo.cxx:193
 TKDNodeInfo.cxx:194
 TKDNodeInfo.cxx:195
 TKDNodeInfo.cxx:196
 TKDNodeInfo.cxx:197
 TKDNodeInfo.cxx:198
 TKDNodeInfo.cxx:199
 TKDNodeInfo.cxx:200
 TKDNodeInfo.cxx:201
 TKDNodeInfo.cxx:202
 TKDNodeInfo.cxx:203
 TKDNodeInfo.cxx:204
 TKDNodeInfo.cxx:205
 TKDNodeInfo.cxx:206
 TKDNodeInfo.cxx:207
 TKDNodeInfo.cxx:208
 TKDNodeInfo.cxx:209
 TKDNodeInfo.cxx:210
 TKDNodeInfo.cxx:211
 TKDNodeInfo.cxx:212
 TKDNodeInfo.cxx:213
 TKDNodeInfo.cxx:214
 TKDNodeInfo.cxx:215
 TKDNodeInfo.cxx:216
 TKDNodeInfo.cxx:217
 TKDNodeInfo.cxx:218
 TKDNodeInfo.cxx:219
 TKDNodeInfo.cxx:220
 TKDNodeInfo.cxx:221
 TKDNodeInfo.cxx:222
 TKDNodeInfo.cxx:223
 TKDNodeInfo.cxx:224
 TKDNodeInfo.cxx:225
 TKDNodeInfo.cxx:226
 TKDNodeInfo.cxx:227
 TKDNodeInfo.cxx:228
 TKDNodeInfo.cxx:229
 TKDNodeInfo.cxx:230
 TKDNodeInfo.cxx:231
 TKDNodeInfo.cxx:232
 TKDNodeInfo.cxx:233
 TKDNodeInfo.cxx:234
 TKDNodeInfo.cxx:235
 TKDNodeInfo.cxx:236
 TKDNodeInfo.cxx:237
 TKDNodeInfo.cxx:238
 TKDNodeInfo.cxx:239
 TKDNodeInfo.cxx:240
 TKDNodeInfo.cxx:241
 TKDNodeInfo.cxx:242
 TKDNodeInfo.cxx:243
 TKDNodeInfo.cxx:244
 TKDNodeInfo.cxx:245
 TKDNodeInfo.cxx:246
 TKDNodeInfo.cxx:247
 TKDNodeInfo.cxx:248
 TKDNodeInfo.cxx:249
 TKDNodeInfo.cxx:250
 TKDNodeInfo.cxx:251
 TKDNodeInfo.cxx:252
 TKDNodeInfo.cxx:253
 TKDNodeInfo.cxx:254
 TKDNodeInfo.cxx:255
 TKDNodeInfo.cxx:256
 TKDNodeInfo.cxx:257
 TKDNodeInfo.cxx:258
 TKDNodeInfo.cxx:259
 TKDNodeInfo.cxx:260
 TKDNodeInfo.cxx:261
 TKDNodeInfo.cxx:262
 TKDNodeInfo.cxx:263
 TKDNodeInfo.cxx:264
 TKDNodeInfo.cxx:265