source: trunk/FACT++/src/DrsCalib.h@ 11678

Last change on this file since 11678 was 11678, checked in by tbretz, 13 years ago
To be unified with the code in MARS.
File size: 7.4 KB
Line 
1class CalibData
2{
3protected:
4 uint64_t fNumEntries;
5 uint64_t fNumCounter;
6
7 size_t fNumSamples;
8 size_t fNumChannels;
9
10 vector<int64_t> fSum;
11 vector<int64_t> fSum2;
12
13 /*
14 vector<map<int16_t, uint8_t> > fMedian;
15 vector<int16_t> fResult;
16 */
17public:
18 CalibData() : fNumEntries(0), fNumSamples(0), fNumChannels(0) { }
19 void Reset()
20 {
21 fNumEntries = 0;
22 fNumSamples = 0;
23 fNumChannels = 0;
24
25 fSum.clear();
26 fSum2.clear();
27 }
28
29 void InitSize(uint16_t channels, uint16_t samples)
30 {
31 fNumChannels = channels;
32 fNumSamples = samples;
33
34 fSum.resize(samples*channels);
35 fSum2.resize(samples*channels);
36
37 //fMedian.resize(40*9*1024);
38 }
39/*
40 continue;
41
42 if (ch<40*9)
43 {
44 fMedian[abs][val[rel]]++;
45
46 int n= 8;
47 if (fNumEntries>0 && fNumEntries%(n*100)==0)
48 {
49 fResult.resize(40*9*1024);
50
51 uint16_t sum = 0;
52 for (map<int16_t, uint8_t>::const_iterator it=fMedian[abs].begin();
53 it!=fMedian[abs].end(); it++)
54 {
55 map<int16_t, uint8_t>::const_iterator is = it;
56 is++;
57 if (is==fMedian[abs].end())
58 break;
59
60 sum += it->second;
61
62 double med = 0;
63
64 med += int64_t(it->second)*int64_t(it->first);
65 med += int64_t(is->second)*int64_t(is->first);
66 med /= int64_t(it->second)+int64_t(is->second);
67
68 if (sum==n*50)
69 {
70 fResult[abs] = it->first;
71 cout << ch << " MED+(50): " << it->first << endl;
72 }
73 if (sum<n*50 && sum+is->second>n*50)
74 {
75 fResult[abs] = med;
76 cout << ch << " MED-(50): " << med << endl;
77 }
78 }
79 cout << ch << " AVG=" << float(fSum[abs])/fNumEntries << endl;
80 fMedian[abs].clear();
81 }
82 } // -2029 -2012 -2003 -1996 -1990 // -1955
83 */
84
85 void AddRel(const int16_t *val, const int16_t *start)
86 {
87 for (size_t ch=0; ch<fNumChannels; ch++)
88 {
89 const int16_t spos = start[ch];
90 if (spos<0)
91 continue;
92
93 const size_t pos = ch*fNumSamples;
94
95 for (size_t i=0; i<fNumSamples; i++)
96 {
97 // Value is relative to trigger
98 // Abs is corresponding index relative to DRS pipeline
99 const size_t rel = pos + i;
100 const size_t abs = pos + (spos+i)%fNumSamples;
101
102 const int64_t v = val[rel];
103
104 fSum[abs] += v;
105 fSum2[abs] += v*v;
106 }
107 }
108
109 fNumEntries++;
110 }
111
112 void AddRel(const int16_t *val, const int16_t *start,
113 const int32_t *offset, const uint32_t scale)
114 {
115 for (size_t ch=0; ch<fNumChannels; ch++)
116 {
117 const int16_t spos = start[ch];
118 if (spos<0)
119 continue;
120
121 const size_t pos = ch*fNumSamples;
122
123 for (size_t i=0; i<fNumSamples; i++)
124 {
125 // Value is relative to trigger
126 // Offset is relative to DRS pipeline
127 // Abs is corresponding index relative to DRS pipeline
128 const size_t rel = pos + i;
129 const size_t abs = pos + (spos+i)%fNumSamples;
130
131 const int64_t v = int64_t(val[rel])*scale-offset[abs];
132
133 fSum[abs] += v;
134 fSum2[abs] += v*v;
135 }
136 }
137
138 fNumEntries++;
139 }
140
141 void AddAbs(const int16_t *val, const int16_t *start,
142 const int32_t *offset, const uint32_t scale)
143 {
144 for (size_t ch=0; ch<fNumChannels; ch++)
145 {
146 const int16_t spos = start[ch];
147 if (spos<0)
148 continue;
149
150 const size_t pos = ch*fNumSamples;
151
152 for (size_t i=0; i<fNumSamples; i++)
153 {
154 // Value is relative to trigger
155 // Offset is relative to DRS pipeline
156 // Abs is corresponding index relative to DRS pipeline
157 const size_t rel = pos + i;
158 const size_t abs = pos + (spos+i)%fNumSamples;
159
160 const int64_t v = int64_t(val[rel])*scale-offset[abs];
161
162 fSum[rel] += v;
163 fSum2[rel] += v*v;
164 }
165 }
166
167 fNumEntries++;
168 }
169
170 static void Apply(int16_t *val, const int16_t *start, uint32_t roi,
171 const int32_t *offset, const uint32_t scaleabs,
172 const int32_t *gain, const uint32_t scalegain,
173 const int32_t *trgoff, const uint32_t scalerel)
174 {
175 for (size_t ch=0; ch<1440; ch++)
176 {
177 const int16_t spos = start[ch];
178 if (spos<0)
179 continue;
180
181 const size_t pos = ch*roi;
182
183 for (size_t i=0; i<roi; i++)
184 {
185 // Value is relative to trigger
186 // Offset is relative to DRS pipeline
187 // Abs is corresponding index relative to DRS pipeline
188 const size_t rel = pos + i;
189 const size_t abs = pos + (spos+i)%1024;
190
191 const int64_t v =
192 + int64_t(val[rel]) *scaleabs*scalerel
193 - int64_t(trgoff[rel])*scaleabs
194 - int64_t(offset[abs])*scalerel
195 ;
196
197 const int64_t div = int64_t(gain[abs])*scaleabs*scalerel;
198
199 val[rel] = (v*scalegain)/div;
200 }
201 }
202 }
203
204 pair<vector<double>,vector<double> > GetSampleStats() const
205 {
206 if (fNumEntries==0)
207 return make_pair(vector<double>(),vector<double>());
208
209 vector<double> mean(fSum.size());
210 vector<double> error(fSum.size());
211
212 vector<int64_t>::const_iterator it = fSum.begin();
213 vector<int64_t>::const_iterator i2 = fSum2.begin();
214 vector<double>::iterator im = mean.begin();
215 vector<double>::iterator ie = error.begin();
216
217 int cnt = 0;
218 while (it!=fSum.end())
219 {
220 *im = /*cnt<fResult.size() ? fResult[cnt] :*/ double(*it)/fNumEntries;
221 *ie = sqrt(double(*i2*int64_t(fNumEntries) - *it * *it))/fNumEntries;
222
223 im++;
224 ie++;
225 it++;
226 i2++;
227 cnt++;
228 }
229
230
231 /*
232 valarray<double> ...
233
234 mean /= fNumEntries;
235 error = sqrt(error/fNumEntries - mean*mean);
236 */
237
238 return make_pair(mean, error);
239 }
240
241 const vector<int64_t> &GetSum() const { return fSum; }
242
243 uint64_t GetNumEntries() const { return fNumEntries; }
244};
245
246
247/*
248Int_t TH1::MyFill(Int_t x, Double_t w)
249{
250 fEntries++;
251
252 fArrayBin[x] += w;
253
254 fSumw2.fArray[x] += w*w
255 fTsumw += w;
256 fTsumw2 += w*w;
257 fTsumwx += w*x;
258 fTsumwx2 += w*x*x;
259 return bin;
260}
261*/
Note: See TracBrowser for help on using the repository browser.