IO factorization for Interpolator
[u/mrichter/AliRoot.git] / STAT / Macros / testkdtreeIF.C
CommitLineData
fe476dd1 1// #include <malloc.h>
2// #include "TMatrixD.h"
3// #include "TRandom.h"
4// #include "TGraph.h"
5// #include "TStopwatch.h"
6//
7
8
9Int_t Mem()
10{
11 // size in 1000 bytes
12 static struct mallinfo memdebug;
13 memdebug = mallinfo();
14 return memdebug.uordblks/1000;
15}
16
17
18
19void testindexes(Int_t nloop, Int_t npoints);
20void testkdtreeIF(Int_t npoints=1000, Int_t nloop=1000, Int_t mode = 1, Int_t bsize=9);
21void TestSizeIF(Int_t npoints=1000, Int_t nrows = 150, Int_t nsec=18, Int_t mode = 1, Int_t bsize=9);
22
23
24
25void testindexes(Int_t nloop, Int_t npoints){
26 //
27 // test indexing
28 //
29 TKDTree<Int_t, Float_t> *kdtree = new TKDTree<Int_t, Float_t>(0,0,0,0);
30 Int_t row =0;
31 Int_t collumn =0;
32 TStopwatch timer;
33 timer.Start();
34 row = 10;
35 for (Int_t iloop=0;iloop<nloop;iloop++)
36 for (Int_t index=1;index<npoints;index++){
37 //row = TMath::Log2(index);
38 //row=0;
39 // if (index< (16<<row)) row=0;
40// for (; index>=(32<<row);row+=5);
41// for (; index>=(2<<row);row++);
42// collumn= index-(1<<row);
43 TKDTree<Int_t, Float_t>::GetCoord(index,row,collumn);
44 //
45 //Int_t index2=kdtree->GetIndex(row,collumn);
46 //printf("%d\t%d\t%d\t%d\n",index,index2,row,collumn);
47 if (kdtree->GetIndex(row,collumn)!=index || collumn<0) {
48 printf("Problem\n");
49 }
50 }
51 timer.Stop();
52 timer.Print();
53}
54
55void TestBuild(Int_t npoints, Int_t bsize){
56 Float_t *data0 = new Float_t[npoints*2];
57 Float_t *data[2];
58 data[0] = &data0[0];
59 data[1] = &data0[npoints];
60 for (Int_t i=0;i<npoints;i++) {
61 data[1][i]= gRandom->Rndm();
62 data[0][i]= gRandom->Rndm();
63 //data[1][i]= 0;
64 //data[0][i]= -i;
65 }
66 Float_t dataf[2];
67 TKDTreeIF *kdtree = new TKDTreeIF(npoints, 2, bsize, data);
68 //kdtree->Build();
69 for (Int_t i=0; i<npoints;i++){
70 Int_t index = kdtree->fIndPoints[i];
71 //printf("%d\t%d\t%f\n",i,index,data[0][index]);
72 }
73 for (Int_t i=0; i<kdtree->fNnodes;i++){
74 //printf("%d\t%f\n",i,kdtree->fNodes[i].fValue);
75 }
76 Float_t sumiter = 0;
77 for (Int_t i=0;i<npoints;i++){
78 dataf[0] = data[0][i];
79 dataf[1] = data[1][i];
80 Int_t index=-1;
81 Int_t iter =0;
82 kdtree->FindPoint(dataf,index, iter);
83 if (i!=index){
84 printf("%d\t%d\t%f\t%f\n",i,index,dataf[0],data[0][index]);
85 }
86 sumiter+=iter;
87 }
88 printf("Mean iter = %f\n",float(sumiter)/float(npoints));
89}
90
91
92void TestSizeIF(Int_t npoints, Int_t nrows, Int_t nsec, Int_t mode, Int_t bsize)
93{
94 //
95 // test size to build kdtree
96 //
97 Int_t before =Mem();
98 for (Int_t isec=0; isec<nsec;isec++)
99 for (Int_t irow=0;irow<nrows;irow++){
100 testkdtreeIF(npoints,0,mode,bsize);
101 }
102 Int_t after = Mem();
103 printf("Memory usage %d\n",after-before);
104}
105
106
107void testkdtreeIF(Int_t npoints, Int_t nloop, Int_t mode, Int_t bsize)
108{
109 //
110 // test speed and functionality of kdtree
111 //
112 Float_t rangey = 100;
113 Float_t rangez = 100;
114 Float_t drangey = 0.1;
115 Float_t drangez = 0.1;
116// Float_t rangey = 20;
117// Float_t rangez = 250;
118// Float_t drangey = 1;
119// Float_t drangez = 1;
120
121 //
122 Float_t *data0 = new Float_t[npoints*2];
123 Float_t *data[2];
124 data[0] = &data0[0];
125 data[1] = &data0[npoints];
126 Int_t i;
127 for (i=0; i<npoints; i++){
128 data[0][i] = gRandom->Uniform(-rangey, rangey);
129 data[1][i] = gRandom->Uniform(-rangez, rangez);
130 // data[i+npoints] = TMath::Nint(gRandom->Uniform(-rangez, rangez)/10.)*10.;
131 //printf("%d %f %f\n", i, data[i], data[i+npoints]);
132 }
133 TStopwatch timerbuild;
134 TKDTree<Int_t, Float_t> *kdtree = new TKDTree<Int_t, Float_t>(npoints,2,bsize,data);
135 kdtree->Build();
136 timerbuild.Stop();
137 timerbuild.Print();
138 TStopwatch timer;
139
140 Float_t countern=0;
141 Float_t counteriter = 0;
142 Float_t counterfound = 0;
143
144
145 if (mode ==2){
146 if (nloop) timer.Start();
147 Int_t res[npoints];
148 Int_t nfound = 0;
149 for (Int_t kloop = 0;kloop<nloop;kloop++){
150 if (kloop==0){
151 counteriter = 0;
152 counterfound= 0;
153 countern = 0;
154 }
155 for (Int_t i=0;i<npoints;i++){
156 Float_t point[2]={data[0][i],data[1][i]};
157 Float_t delta[2]={drangey,drangez};
158 Int_t iter =0;
159 nfound =0;
160 Int_t bnode =0;
161 //kdtree->FindBNode(point,delta, bnode);
162 //continue;
163 kdtree->FindInRangeA(point,delta,res,nfound,iter,bnode);
164 if (kloop==0){
165 //Bool_t isOK = kTRUE;
166 Bool_t isOK = kFALSE;
167 for (Int_t ipoint=0;ipoint<nfound;ipoint++)
168 if (res[ipoint]==i) isOK =kTRUE;
169 counteriter+=iter;
170 counterfound+=nfound;
171 if (isOK) {
172 countern++;
173 }else{
174 printf("Bug\n");
175 }
176 }
177 }
178 }
179 if (nloop){
180 timer.Stop();
181 timer.Print();
182 }
183 }
184 delete [] data0;
185
186 counteriter/=npoints;
187 counterfound/=npoints;
188 if (nloop) printf("Find nearest point:\t%f\t%f\t%f\n",countern, counteriter, counterfound);
189}