1 | #include "Sample.h"
|
---|
2 |
|
---|
3 | Sample::Sample()
|
---|
4 | {
|
---|
5 | InitMembers();
|
---|
6 | mpRndNumList = new vector<int>;
|
---|
7 | }
|
---|
8 |
|
---|
9 | Sample::Sample(int size)
|
---|
10 | {
|
---|
11 | InitMembers();
|
---|
12 | mpRndNumList = new vector<int>;
|
---|
13 | mpRndNumList->resize(size);
|
---|
14 | }
|
---|
15 |
|
---|
16 | Sample::Sample(vector<int> *RndNumList)
|
---|
17 | {
|
---|
18 | InitMembers();
|
---|
19 | mpRndNumList = RndNumList;
|
---|
20 | }
|
---|
21 |
|
---|
22 | Sample::Sample(vector<int> *RndNumList, int size)
|
---|
23 | {
|
---|
24 | InitMembers();
|
---|
25 | mpRndNumList = RndNumList;
|
---|
26 | mpRndNumList->resize(size);
|
---|
27 | }
|
---|
28 |
|
---|
29 | Sample::~Sample(void)
|
---|
30 | {
|
---|
31 | if (mpRndNumList != NULL)
|
---|
32 | {
|
---|
33 | delete mpRndNumList;
|
---|
34 | }
|
---|
35 | }
|
---|
36 |
|
---|
37 | void Sample::InitMembers()
|
---|
38 | {
|
---|
39 | mpRndNumList = NULL;
|
---|
40 | mMinNumber = 0;
|
---|
41 | mMaxNumber = 0;
|
---|
42 | mSampleSize = 0;
|
---|
43 | mVerbosityLvl = 0;
|
---|
44 | mSeed = 4357;
|
---|
45 | }
|
---|
46 |
|
---|
47 | //============================================================================
|
---|
48 | //ACCESS
|
---|
49 | //============================================================================
|
---|
50 | void Sample::SetMinNumber(int min){ mMinNumber = min; return;}
|
---|
51 | void Sample::SetMaxNumber(int max){ mMaxNumber = max; return;}
|
---|
52 | void Sample::SetSampleSize(int size){mSampleSize = size; return;}
|
---|
53 | void Sample::SetSeed(int seed)
|
---|
54 | {
|
---|
55 | mSeed = seed;
|
---|
56 | mRandom.SetSeed(mSeed);
|
---|
57 | return;
|
---|
58 | }
|
---|
59 | void Sample::SetVerbosityLvl(int VerbosityLvl){ mVerbosityLvl = VerbosityLvl; return;}
|
---|
60 | int Sample::GetMinNumber(){return mMinNumber;}
|
---|
61 | int Sample::GetMaxNumber(){return mMaxNumber;}
|
---|
62 | int Sample::GetSampleSize(){return mSampleSize;}
|
---|
63 | int Sample::GetSeed(){return mSeed;}
|
---|
64 |
|
---|
65 | void Sample::GetRndListValues(vector<int> *RndNumList)
|
---|
66 | {
|
---|
67 | vector<int>::iterator it;
|
---|
68 |
|
---|
69 | for ( it=mpRndNumList->begin() ; it < mpRndNumList->end(); it++ )
|
---|
70 | {
|
---|
71 | RndNumList->push_back(*it);
|
---|
72 | }
|
---|
73 |
|
---|
74 | return;
|
---|
75 | }
|
---|
76 |
|
---|
77 | vector<int>* Sample::GetRndListPtr()
|
---|
78 | {
|
---|
79 | return mpRndNumList;
|
---|
80 | }
|
---|
81 |
|
---|
82 | //============================================================================
|
---|
83 | //OPERATIONS
|
---|
84 | //============================================================================
|
---|
85 | int
|
---|
86 | Sample::GenerateRandomInt(
|
---|
87 | int min,
|
---|
88 | int max)
|
---|
89 | {
|
---|
90 | // int rndNum = min + mRandom.Integer(max - min);
|
---|
91 | int rndNum = (int) mRandom.Uniform(min, max);
|
---|
92 | // cout << rndNum << endl;
|
---|
93 | return rndNum;
|
---|
94 | }
|
---|
95 |
|
---|
96 | int
|
---|
97 | Sample::GenerateRandomInt()
|
---|
98 | {
|
---|
99 | return GenerateRandomInt( mMinNumber, mMaxNumber );
|
---|
100 | // return 0;
|
---|
101 | }
|
---|
102 |
|
---|
103 | void
|
---|
104 | Sample::GenerateSample()
|
---|
105 | {
|
---|
106 | GenerateSample( mpRndNumList, mMinNumber, mMaxNumber, mSampleSize);
|
---|
107 | return;
|
---|
108 | }
|
---|
109 |
|
---|
110 | void
|
---|
111 | Sample::GenerateSample(
|
---|
112 | vector<int>* sampleVector,
|
---|
113 | int min,
|
---|
114 | int max,
|
---|
115 | int size)
|
---|
116 | {
|
---|
117 | if ( size > max - min)
|
---|
118 | {
|
---|
119 | cout << "sample size is larger than range of sample, will set size to range" << endl;
|
---|
120 | size = max - min + 1;
|
---|
121 | }
|
---|
122 | //resize destination vector for pulled numbers
|
---|
123 | sampleVector->resize(size);
|
---|
124 |
|
---|
125 | //calculate qunatity of numbers in range
|
---|
126 | int qunatityOfNumbers = max - min + 1;
|
---|
127 |
|
---|
128 | //crate a vector list of ordered numbers in defined range
|
---|
129 | vector<int> listOfNumbers;
|
---|
130 | listOfNumbers.resize(qunatityOfNumbers);
|
---|
131 |
|
---|
132 | //fill a list of ordered numbers in defined range
|
---|
133 | for (int i = min; i <= max; i++)
|
---|
134 | {
|
---|
135 | listOfNumbers.at(i)=i;
|
---|
136 | }
|
---|
137 |
|
---|
138 | //container to fill in numbers with ordering
|
---|
139 | set<int> sampleSet;
|
---|
140 |
|
---|
141 | //container for insert result
|
---|
142 | bool result;
|
---|
143 | for (int i = 0; i < size; i++)
|
---|
144 | {
|
---|
145 | int randomNumber = GenerateRandomInt( 0, qunatityOfNumbers-i);
|
---|
146 | result = sampleSet.insert(listOfNumbers.at(randomNumber)).second;
|
---|
147 | if (result)
|
---|
148 | {
|
---|
149 | // cout << "rndNR " << randomNumber << endl;
|
---|
150 | listOfNumbers.erase(listOfNumbers.begin()+randomNumber);
|
---|
151 | }
|
---|
152 | else if (!result)
|
---|
153 | {
|
---|
154 | cout << " pulled number exists, pulling again" << endl;
|
---|
155 | i--;
|
---|
156 | }
|
---|
157 | }
|
---|
158 |
|
---|
159 | set<int>::iterator it;
|
---|
160 |
|
---|
161 | int counter = 0;
|
---|
162 | for ( it=sampleSet.begin() ; it != sampleSet.end(); it++ )
|
---|
163 | {
|
---|
164 | sampleVector->at(counter)=*it;
|
---|
165 | counter++;
|
---|
166 | }
|
---|
167 |
|
---|
168 | return;
|
---|
169 | }
|
---|
170 |
|
---|
171 |
|
---|
172 | int
|
---|
173 | Sample::BootstrapVector(vector<double> *inVector, vector<double> *outVector)
|
---|
174 | {
|
---|
175 | //get size of sample to be bootstrapped
|
---|
176 | int sampleSize = inVector->size();
|
---|
177 |
|
---|
178 | //Vector with positions in original sample
|
---|
179 | vector<int> entryID (sampleSize,0);
|
---|
180 |
|
---|
181 | //calculate wich entries from inVector will be put into outVector
|
---|
182 | BootstrapSample(&entryID, 0, sampleSize, sampleSize );
|
---|
183 |
|
---|
184 | vector<int>::iterator it;
|
---|
185 |
|
---|
186 | //Loop over entryID vector to fill content from inVector to outVector
|
---|
187 | int counter = 0;
|
---|
188 | for ( it=entryID.begin() ; it != entryID.end(); it++ )
|
---|
189 | {
|
---|
190 | outVector->at(counter) = inVector->at(*it);
|
---|
191 | counter++;
|
---|
192 | }
|
---|
193 |
|
---|
194 | return counter + 1;
|
---|
195 | }
|
---|
196 |
|
---|
197 | int
|
---|
198 | Sample::BootstrapTH1(TH1* inputHisto, TH1* outHisto)
|
---|
199 | {
|
---|
200 | //reset outHisto in case it is reused
|
---|
201 | outHisto->Reset();
|
---|
202 |
|
---|
203 | //compute number of bins for 1-d histogram h1
|
---|
204 | int nbins = inputHisto->GetXaxis()->GetNbins();
|
---|
205 |
|
---|
206 | if (mVerbosityLvl > 1) cout << "nbins: " << nbins << endl;
|
---|
207 | //we need to get the binning
|
---|
208 |
|
---|
209 | //vector to contain entries of TH1
|
---|
210 | vector<double> entries;
|
---|
211 |
|
---|
212 | //quantity of entries in bin
|
---|
213 | int quantity = 0; //number of entries in a bin
|
---|
214 | double value = 0; //value of a bin
|
---|
215 |
|
---|
216 | //Loop over bins to fill entries vektor
|
---|
217 | for (int i=1;i<nbins;i++)
|
---|
218 | {
|
---|
219 | //value = inputHisto->GetBinLowEdge(i);
|
---|
220 | value = inputHisto->GetBinCenter(i);
|
---|
221 | quantity = inputHisto->GetBinContent(i);
|
---|
222 | //Loop over bin quantities
|
---|
223 | for (int j = 0; j < quantity; j++)
|
---|
224 | {
|
---|
225 | // fill entries vektor with value of a bin
|
---|
226 | // as many times as the bin has entries
|
---|
227 | entries.push_back(value);
|
---|
228 | }
|
---|
229 | }
|
---|
230 |
|
---|
231 | //get size of sample to be bootstrapped
|
---|
232 | int sampleSize = entries.size();
|
---|
233 |
|
---|
234 | if (mVerbosityLvl > 1) cout << "sampleSize: " << sampleSize << endl;
|
---|
235 |
|
---|
236 | //Vector with positions in original sample
|
---|
237 | vector<int> entryID (sampleSize,0);
|
---|
238 |
|
---|
239 | //calculate a list with random int numbers between 0 ad sampleSize
|
---|
240 | // and fill it into entryID vector
|
---|
241 | BootstrapSample(&entryID, 0, sampleSize, sampleSize );
|
---|
242 |
|
---|
243 | //Loop over entryID vector to bootstrap the histogram
|
---|
244 | int counter = 0;
|
---|
245 | for ( unsigned int i = 0 ; i < entryID.size(); i++ )
|
---|
246 | {
|
---|
247 | //fill values pulled from entries vektor into outHisto
|
---|
248 | outHisto->Fill(
|
---|
249 | entries.at( entryID.at(i) )
|
---|
250 | );
|
---|
251 | counter++;
|
---|
252 | }
|
---|
253 |
|
---|
254 | //return the nummber of filled entries
|
---|
255 | return counter + 1;
|
---|
256 | }
|
---|
257 |
|
---|
258 | void
|
---|
259 | Sample::BootstrapSample()
|
---|
260 | {
|
---|
261 | BootstrapSample( mpRndNumList, mMinNumber, mMaxNumber, mSampleSize );
|
---|
262 |
|
---|
263 | return;
|
---|
264 | }
|
---|
265 |
|
---|
266 | void
|
---|
267 | Sample::BootstrapSample(
|
---|
268 | vector<int> *sampleVector,
|
---|
269 | int numMinEvent,
|
---|
270 | int numMaxEvent,
|
---|
271 | int size)
|
---|
272 | {
|
---|
273 | if (size == 0){
|
---|
274 | if (mVerbosityLvl > 1) cout << "Bootstrapping: size of vector = 0; nothing to do" << endl;
|
---|
275 | return;
|
---|
276 | }
|
---|
277 |
|
---|
278 | //resize the sample vector to size of boostrapped sample
|
---|
279 | sampleVector->resize(size);
|
---|
280 |
|
---|
281 | //list of rndnumbers pulled with putting back
|
---|
282 | multiset<int> sampleSet;
|
---|
283 |
|
---|
284 | //loop over samplesize to generate random numbers to fill into sampleset
|
---|
285 | for (int i = 0; i < size; i++)
|
---|
286 | {
|
---|
287 | int randomNumber = GenerateRandomInt( numMinEvent, numMaxEvent);
|
---|
288 | sampleSet.insert(randomNumber);
|
---|
289 | }
|
---|
290 |
|
---|
291 | multiset<int>::iterator it;
|
---|
292 |
|
---|
293 | // set<int>::iterator it;
|
---|
294 |
|
---|
295 | int counter = 0;
|
---|
296 | // loop over list of rndnumbers and fill their entries into sampleVector
|
---|
297 | // entries are vector positions
|
---|
298 | for ( it=sampleSet.begin() ; it != sampleSet.end(); it++ )
|
---|
299 | {
|
---|
300 | sampleVector->at(counter)=*it;
|
---|
301 | counter++;
|
---|
302 | }
|
---|
303 |
|
---|
304 | return;
|
---|
305 | }
|
---|
306 | //NOTE: crashes for some reason if size is smaller than max-min
|
---|