source: trunk/FACT++/src/feedback.cc@ 17203

Last change on this file since 17203 was 17203, checked in by tbretz, 11 years ago
Added custom channel offsets to correct gain variations seen in the camera; some small imporvements to the output; instead of using Ubd as reference, now the voltage around which the gain is spposed to be 0 is used as reference - this also improves the behaviour during bright light.
File size: 41.9 KB
Line 
1#include <valarray>
2#include <algorithm>
3
4#include "Dim.h"
5#include "Event.h"
6#include "Shell.h"
7#include "StateMachineDim.h"
8#include "Connection.h"
9#include "Configuration.h"
10#include "Console.h"
11#include "externals/PixelMap.h"
12#include "externals/Interpolator2D.h"
13
14#include "tools.h"
15
16#include "LocalControl.h"
17
18#include "HeadersFSC.h"
19#include "HeadersBIAS.h"
20#include "HeadersFeedback.h"
21
22#include "DimState.h"
23#include "DimDescriptionService.h"
24
25using namespace std;
26
27// ------------------------------------------------------------------------
28
29class StateMachineFeedback : public StateMachineDim
30{
31private:
32 PixelMap fMap;
33
34 bool fIsVerbose;
35
36 DimVersion fDim;
37
38 DimDescribedState fDimFSC;
39 DimDescribedState fDimBias;
40
41 DimDescribedService fDimCalibration;
42 DimDescribedService fDimCalibration2;
43 DimDescribedService fDimCalibrationR8;
44 DimDescribedService fDimCurrents;
45 DimDescribedService fDimOffsets;
46
47 vector<float> fCalibCurrentMes[6]; // Measured calibration current at six different levels
48 vector<float> fCalibVoltage[6]; // Corresponding voltage as reported by biasctrl
49
50 vector<int64_t> fCurrentsAvg;
51 vector<int64_t> fCurrentsRms;
52
53 vector<float> fVoltGapd; // Nominal breakdown voltage + 1.1V
54 vector<float> fBiasVolt; // Output voltage as reported by bias crate (voltage between R10 and R8)
55 vector<float> fBiasR9; //
56 vector<uint16_t> fBiasDac; // Dac value corresponding to the voltage setting
57
58 vector<float> fCalibration;
59 vector<float> fCalibDeltaI;
60 vector<float> fCalibR8;
61
62 int64_t fCursorCur;
63
64 Time fTimeCalib;
65 Time fTimeTemp;
66
67 double fUserOffset;
68 vector<double> fTempOffset;
69 float fTempOffsetAvg;
70 float fTempOffsetRms;
71 double fTempCoefficient;
72 double fTemp;
73
74 vector<double> fVoltOffset;
75
76 uint16_t fCurrentRequestInterval;
77 uint16_t fNumCalibIgnore;
78 uint16_t fNumCalibRequests;
79 uint16_t fCalibStep;
80
81 // ============================= Handle Services ========================
82
83 int HandleBiasStateChange()
84 {
85 if (fDimBias.state()==BIAS::State::kVoltageOn && GetCurrentState()==Feedback::State::kCalibrating)
86 {
87 Dim::SendCommandNB("BIAS_CONTROL/REQUEST_STATUS");
88 Info("Starting calibration step "+to_string(fCalibStep));
89 }
90
91 if (fDimBias.state()==BIAS::State::kVoltageOff && GetCurrentState()==Feedback::State::kInProgress)
92 return Feedback::State::kCalibrated;
93
94 return GetCurrentState();
95 }
96 // ============================= Handle Services ========================
97
98 bool CheckEventSize(size_t has, const char *name, size_t size)
99 {
100 if (has==size)
101 return true;
102
103 // Disconnected
104 if (has==0)
105 return false;
106
107 ostringstream msg;
108 msg << name << " - Received event has " << has << " bytes, but expected " << size << ".";
109 Fatal(msg);
110 return false;
111 }
112
113 int HandleBiasNom(const EventImp &evt)
114 {
115 if (evt.GetSize()>=416*sizeof(float))
116 {
117 fVoltGapd.assign(evt.Ptr<float>(), evt.Ptr<float>()+416);
118 fBiasR9.assign(evt.Ptr<float>()+2*416, evt.Ptr<float>()+3*416);
119 Info("Nominal bias voltages and calibration resistor received.");
120 }
121
122 return GetCurrentState();
123 }
124
125 int HandleBiasVoltage(const EventImp &evt)
126 {
127 if (evt.GetSize()>=416*sizeof(float))
128 fBiasVolt.assign(evt.Ptr<float>(), evt.Ptr<float>()+416);
129 return GetCurrentState();
130 }
131
132 int HandleBiasDac(const EventImp &evt)
133 {
134 if (evt.GetSize()>=416*sizeof(uint16_t))
135 fBiasDac.assign(evt.Ptr<uint16_t>(), evt.Ptr<uint16_t>()+416);
136 return GetCurrentState();
137 }
138
139 int HandleCameraTemp(const EventImp &evt)
140 {
141 if (!CheckEventSize(evt.GetSize(), "HandleCameraTemp", 323*sizeof(float)))
142 {
143 fTimeTemp = Time(Time::none);
144 return GetCurrentState();
145 }
146
147 //fTempOffset = (avgt-25)*0.0561765; // [V] From Hamamatsu datasheet
148 //fTempOffset = (avgt-25)*0.05678; // [V] From Hamamatsu datasheet plus our own measurement (gein vs. temperature)
149
150 const float *ptr = evt.Ptr<float>(4);
151
152 fTimeTemp = evt.GetTime();
153 fTemp = evt.Get<float>(321*4);
154
155 fTempOffsetAvg = (fTemp-25)*fTempCoefficient;
156 fTempOffsetRms = evt.Get<float>(322*4)*fTempCoefficient;
157
158 fTempOffset.resize(320);
159 for (int i=0; i<320; i++)
160 fTempOffset[i] = (ptr[i]-25)*fTempCoefficient;
161
162 return GetCurrentState();
163 }
164
165 pair<vector<float>, vector<float>> AverageCurrents(const int16_t *ptr, int n)
166 {
167 if (fCursorCur++>=0)
168 {
169 for (int i=0; i<BIAS::kNumChannels; i++)
170 {
171 fCurrentsAvg[i] += ptr[i];
172 fCurrentsRms[i] += ptr[i]*ptr[i];
173 }
174 }
175
176 if (fCursorCur<n)
177 return make_pair(vector<float>(), vector<float>());
178
179 const double conv = 5e-3/4096;
180
181 vector<float> rms(BIAS::kNumChannels);
182 vector<float> avg(BIAS::kNumChannels);
183 for (int i=0; i<BIAS::kNumChannels; i++)
184 {
185 avg[i] = double(fCurrentsAvg[i])/fCursorCur * conv;
186 rms[i] = double(fCurrentsRms[i])/fCursorCur * conv * conv;
187 rms[i] -= avg[i]*avg[i];
188 rms[i] = rms[i]<0 ? 0 : sqrt(rms[i]);
189 }
190
191 return make_pair(avg, rms);
192 }
193
194 int HandleCalibration(const EventImp &evt)
195 {
196 if (fDimBias.state()!=BIAS::State::kVoltageOn)
197 return GetCurrentState();
198
199 const uint16_t dac = 256+512*fCalibStep; // Command value
200
201 // Only the channels which are no spare channels are ramped
202 if (std::count(fBiasDac.begin(), fBiasDac.end(), dac)!=320)
203 return GetCurrentState();
204
205 const auto rc = AverageCurrents(evt.Ptr<int16_t>(), fNumCalibRequests);
206 if (rc.first.size()==0)
207 {
208 Dim::SendCommandNB("BIAS_CONTROL/REQUEST_STATUS");
209 return GetCurrentState();
210 }
211
212 const vector<float> &avg = rc.first;
213 const vector<float> &rms = rc.second;
214
215 // Current through resistor R8
216 fCalibCurrentMes[fCalibStep] = avg; // [A]
217 fCalibVoltage[fCalibStep] = fBiasVolt; // [V]
218
219 // ------------------------- Update calibration data --------------------
220
221 struct cal_data
222 {
223 uint32_t dac;
224 float U[416];
225 float Iavg[416];
226 float Irms[416];
227
228 cal_data() { memset(this, 0, sizeof(cal_data)); }
229 } __attribute__((__packed__));
230
231 cal_data cal;
232 cal.dac = dac;
233 memcpy(cal.U, fBiasVolt.data(), 416*sizeof(float));
234 memcpy(cal.Iavg, avg.data(), 416*sizeof(float));
235 memcpy(cal.Irms, rms.data(), 416*sizeof(float));
236
237 fDimCalibration2.setData(cal);
238 fDimCalibration2.Update(fTimeCalib);
239
240 // -------------------- Start next calibration steo ---------------------
241
242 if (++fCalibStep<6)
243 {
244 fCursorCur = -fNumCalibIgnore;
245 fCurrentsAvg.assign(BIAS::kNumChannels, 0);
246 fCurrentsRms.assign(BIAS::kNumChannels, 0);
247
248 Dim::SendCommandNB("BIAS_CONTROL/SET_GLOBAL_DAC", uint32_t(256+512*fCalibStep));
249
250 return GetCurrentState();
251 }
252
253 // --------------- Calculate old style calibration ----------------------
254
255 fCalibration.resize(BIAS::kNumChannels*4);
256
257 float *pavg = fCalibration.data();
258 float *prms = fCalibration.data()+BIAS::kNumChannels;
259 float *pres = fCalibration.data()+BIAS::kNumChannels*2;
260 float *pUmes = fCalibration.data()+BIAS::kNumChannels*3;
261
262 for (int i=0; i<BIAS::kNumChannels; i++)
263 {
264 const double I = fCalibCurrentMes[5][i]; // [A]
265 const double U = fBiasVolt[i]; // [V]
266
267 pavg[i] = I*1e6; // [uA]
268 prms[i] = rms[i]*1e6; // [uA]
269 pres[i] = U/I; // [Ohm]
270 pUmes[i] = U; // [V]
271 }
272
273 fDimCalibration.setData(fCalibration);
274 fDimCalibration.Update(fTimeCalib);
275
276 // -------------------- New style calibration --------------------------
277
278 fCalibDeltaI.resize(BIAS::kNumChannels);
279 fCalibR8.resize(BIAS::kNumChannels);
280
281 // Linear regression of the values at 256+512*N for N={ 3, 4, 5 }
282 for (int i=0; i<BIAS::kNumChannels; i++)
283 {
284 // x: Idac
285 // y: Iadc
286
287 double x = 0;
288 double y = 0;
289 double xx = 0;
290 double xy = 0;
291
292 const int beg = 3;
293 const int end = 5;
294 const int len = end-beg+1;
295
296 for (int j=beg; j<=end; j++)
297 {
298 const double Idac = (256+512*j)*1e-3/4096;
299
300 x += Idac;
301 xx += Idac*Idac;
302 y += fCalibCurrentMes[j][i];
303 xy += fCalibCurrentMes[j][i]*Idac;
304 }
305
306 const double m1 = xy - x*y / len;
307 const double m2 = xx - x*x / len;
308
309 const double m = m2==0 ? 0 : m1/m2;
310
311 const double t = (y - m*x) / len;
312
313 fCalibDeltaI[i] = t; // [A]
314 fCalibR8[i] = 100/m; // [Ohm]
315 }
316
317 vector<float> v;
318 v.reserve(BIAS::kNumChannels*2);
319 v.insert(v.end(), fCalibDeltaI.begin(), fCalibDeltaI.end());
320 v.insert(v.end(), fCalibR8.begin(), fCalibR8.end());
321
322 fDimCalibrationR8.setData(v);
323 fDimCalibrationR8.Update(fTimeCalib);
324
325 // ---------------------------------------------------------------------
326
327 Info("Calibration successfully done.");
328 Dim::SendCommandNB("BIAS_CONTROL/SET_ZERO_VOLTAGE");
329
330 return Feedback::State::kCalibrated;
331 }
332
333 int HandleBiasCurrent(const EventImp &evt)
334 {
335 if (!CheckEventSize(evt.GetSize(), "HandleBiasCurrent", BIAS::kNumChannels*sizeof(uint16_t)))
336 return Feedback::State::kConnected;
337
338 if (GetCurrentState()<Feedback::State::kCalibrating)
339 return GetCurrentState();
340
341 // ------------------------------- HandleCalibration -----------------------------------
342 if (GetCurrentState()==Feedback::State::kCalibrating)
343 return HandleCalibration(evt);
344
345 // ---------------------- Calibrated, WaitingForData, InProgress -----------------------
346
347 // We are waiting but no valid temperature yet, go on waiting
348 if (GetCurrentState()==Feedback::State::kWaitingForData &&
349 (!fTimeTemp.IsValid() || Time()-fTimeTemp>boost::posix_time::minutes(5)))
350 return GetCurrentState();
351
352 // We are already in progress but no valid temperature update anymore
353 if (GetCurrentState()==Feedback::State::kInProgress &&
354 (!fTimeTemp.IsValid() || Time()-fTimeTemp>boost::posix_time::minutes(5)))
355 {
356 Warn("Current control in progress, but last received temperature older than 5min... switching voltage off.");
357 Dim::SendCommandNB("BIAS_CONTROL/SET_ZERO_VOLTAGE");
358 return Feedback::State::kCalibrated;
359 }
360
361 // ---------------------- Calibrated, WaitingForData, InProgress -----------------------
362
363 const int Navg = fDimBias.state()!=BIAS::State::kVoltageOn ? 1 : 3;
364
365 const vector<float> &Imes = AverageCurrents(evt.Ptr<int16_t>(), Navg).first;
366 if (Imes.size()==0)
367 return GetCurrentState();
368
369 fCurrentsAvg.assign(416, 0);
370 fCurrentsRms.assign(416, 0);
371 fCursorCur = 0;
372
373 // -------------------------------------------------------------------------------------
374
375 // Nominal overvoltage (w.r.t. the bias setup values)
376 const double overvoltage = GetCurrentState()<Feedback::State::kWaitingForData ? 0 : fUserOffset;
377
378 double avg[2] = { 0, 0 };
379 double min[2] = { 90, 90 };
380 double max[2] = { -90, -90 };
381 int num[3] = { 0, 0, 0 };
382
383 vector<double> med[3];
384 med[0].resize(416);
385 med[1].resize(416);
386 med[2].resize(416);
387
388 struct dim_data
389 {
390 float I[416];
391 float Iavg;
392 float Irms;
393 float Imed;
394 float Idev;
395 uint32_t N;
396 float Tdiff;
397 float Uov[416];
398 float Unom;
399 float dUtemp;
400
401 dim_data() { memset(this, 0, sizeof(dim_data)); }
402 } __attribute__((__packed__));
403
404 int Ndev[3] = { 0, 0, 0 };
405
406 dim_data data;
407
408 data.Unom = overvoltage;
409 data.dUtemp = fTempOffsetAvg;
410
411 vector<float> vec(416);
412
413 /*
414 if (fEnableOldAlgorithm)
415 {
416 // ================================= old =======================
417 // Pixel 583: 5 31 == 191 (5) C2 B3 P3
418 // Pixel 830: 2 2 == 66 (4) C0 B8 P1
419 // Pixel 1401: 6 1 == 193 (5) C2 B4 P0
420
421 // 3900 Ohm/n + 1000 Ohm + 1100 Ohm (with n=4 or n=5)
422 const double R[2] = { 3075, 2870 };
423
424 const float *Iavg = fCalibration.data(); // Offset at U=fCalibrationOffset
425 const float *Ravg = fCalibration.data()+BIAS::kNumChannels*2; // Measured resistance
426
427 for (int i=0; i<320; i++)
428 {
429 const PixelMapEntry &hv = fMap.hv(i);
430 if (!hv)
431 continue;
432
433 // Average measured current
434 const double Im = Imes[i] * 1e6; // [uA]
435
436 // Group index (0 or 1) of the of the pixel (4 or 5 pixel patch)
437 const int g = hv.group();
438
439 // Serial resistors in front of the G-APD
440 double Rg = R[g];
441
442 // This is assuming that the broken pixels have a 390 Ohm instead of 3900 Ohm serial resistor
443 if (i==66) // Pixel 830(66)
444 Rg = 2400; // 2400 = (3/3900 + 1/390) + 1000 + 1100
445 if (i==191 || i==193) // Pixel 583(191) / Pixel 1401(193)
446 Rg = 2379; // 2379 = (4/3900 + 1/390) + 1000 + 1100
447
448 const double r = 1./(1./Rg + 1./Ravg[i]); // [Ohm]
449
450 // Offset induced by the voltage above the calibration point
451 const double Ubd = fVoltGapd[i] + fTempOffsets[i];
452 const double U0 = Ubd + overvoltage - fCalibVoltage[5][i]; // appliedOffset-fCalibrationOffset;
453 const double dI = U0/Ravg[i]; // [V/Ohm]
454
455 // Offset at the calibration point (make sure that the calibration is
456 // valid (Im[i]>Iavg[i]) and we operate above the calibration point)
457 const double I = Im>Iavg[i] ? Im - Iavg[i] : 0; // [uA]
458
459 // Make sure that the averaged resistor is valid
460 const double dU = Ravg[i]>10000 ? r*(I*1e-6 - dI) : 0;
461
462 vec[i] = Ubd + overvoltage + dU;
463
464 // Calculate statistics only for channels with a valid calibration
465 if (Iavg[i]>0)
466 {
467 med[g][num[g]] = dU;
468 avg[g] += dU;
469 num[g]++;
470
471 if (dU<min[g])
472 min[g] = dU;
473 if (dU>max[g])
474 max[g] = dU;
475
476 data.I[i] = Imes[i]*1e6 - fBiasVolt[i]/Ravg[i]*1e6;
477 data.I[i] /= hv.count();
478
479 if (i==66)
480 data.I[i] /= 1.3;
481 if (i==191 || i==193)
482 data.I[i] /= 1.4;
483
484 data.Iavg += data.I[i];
485 data.Irms += data.I[i]*data.I[i];
486
487 med[2][num[2]++] = data.I[i];
488 }
489 }
490 }
491 */
492
493 for (int i=0; i<320/*BIAS::kNumChannels*/; i++)
494 {
495 const PixelMapEntry &hv = fMap.hv(i);
496 if (!hv)
497 continue;
498
499 // Number of G-APDs in this patch
500 const int N = hv.count();
501
502 // Average measured ADC value for this channel
503 const double adc = Imes[i]/* * (5e-3/4096)*/; // [A]
504
505 // Current through ~100 Ohm measurement resistor
506 const double I8 = (adc-fCalibDeltaI[i])*fCalibR8[i]/100;
507
508 // Current through calibration resistors (R9)
509 // This is uncalibrated, biut since the corresponding calibrated
510 // value I8 is subtracted, the difference should yield a correct value
511 const double I9 = fBiasDac[i] * (1e-3/4096);//U9/R9; [A]
512
513 // Current in R4/R5 branch
514 const double Iout = I8 - I9;//I8>I9 ? I8 - I9 : 0;
515
516 // Applied voltage at calibration resistors, according to biasctrl
517 const double U9 = fBiasVolt[i];
518
519 // Serial resistors (one 1kOhm at the output of the bias crate, one 1kOhm in the camera)
520 const double R4 = 2000;
521
522 // Serial resistor of the individual G-APDs
523 double R5 = 3900./N;
524
525 // This is assuming that the broken pixels have a 390 Ohm instead of 3900 Ohm serial resistor
526 if (i==66) // Pixel 830(66)
527 R5 = 300; // 2400 = 1/(3/3900 + 1/390)
528 if (i==191 || i==193) // Pixel 583(191) / Pixel 1401(193)
529 R5 = 390/1.4; // 379 = 1/(4/3900 + 1/390)
530
531 // The measurement resistor
532 const double R8 = 100;
533
534 // Total resistance of branch with diodes (R4+R5)
535 // Assuming that the voltage output of the OpAMP is linear
536 // with the DAC setting and not the voltage at R9, the
537 // additional voltage drop at R8 must be taken into account
538 const double R = R4 + R5 + R8;
539
540 // For the patches with a broken resistor - ignoring the G-APD resistance -
541 // we get:
542 //
543 // I[R=3900] = Iout * 1/(10+(N-1)) = Iout /(N+9)
544 // I[R= 390] = Iout * (1 - 1/ (10+(N-1))) = Iout * (N+8)/(N+9)
545 //
546 // I[R=390] / I[R=3900] = N+8
547 //
548 // Udrp = Iout*3900/(N+9) + Iout*1000 + Iout*1000 = Iout * R
549
550 // Voltage drop in R4/R5 branch (for the G-APDs with correct resistor)
551 const double Udrp = R*Iout;
552
553 // Nominal breakdown voltage with correction for temperature dependence
554 const double Ubd = fVoltGapd[i] + fVoltOffset[i] + fTempOffset[i];
555
556 // Current overvoltage (at a G-APD with the correct 3900 Ohm resistor)
557 //const double Uov = U9-Udrp-Ubd>0 ? U9-Udrp-Ubd : 0;
558 const double Uov = U9-Udrp-Ubd>-0.34 ? U9-Udrp-Ubd : -0.34;
559
560 // Iout linear with U9 above Ubd
561 //
562 // Rx = (U9-Ubd)/Iout
563 // I' = (U9'-Ubd) / Rx
564 // Udrp' = R*I'
565 // Uov = U9' - Udrp' - Ubd
566 // Uov = overvoltage
567 //
568 // overvoltage = U9' - Udrp' - Ubd
569 // overvoltage = U9' - R*I' - Ubd
570 // overvoltage = U9' - R*((U9'-Ubd)/Rx) - Ubd
571 // overvoltage = U9' - U9'*R/Rx + Ubd*R/Rx - Ubd
572 // overvoltage = U9'*(1 - R/Rx) + Ubd*R/Rx - Ubd
573 // overvoltage - Ubd*R/Rx +Ubd = U9'*(1 - R/Rx)
574 // U9' = [ overvoltage - Ubd*R/Rx +Ubd ] / (1 - R/Rx)
575 //
576
577 // The current through one G-APD is the sum divided by the number of G-APDs
578 // (assuming identical serial resistors)
579 double Iapd = Iout/N;
580
581 // In this and the previosu case we neglect the resistance of the G-APDs, but we can make an
582 // assumption: The differential resistance depends more on the NSB than on the PDE,
583 // thus it is at least comparable for all G-APDs in the patch. In addition, although the
584 // G-APD with the 390Ohm serial resistor has the wrong voltage applied, this does not
585 // significantly influences the ohmic resistor or the G-APD because the differential
586 // resistor is large enough that the increase of the overvoltage does not dramatically
587 // increase the current flow as compared to the total current flow.
588 if (i==66 || i==191 || i==193)
589 Iapd = Iout/(N+9); // Iapd = R5*Iout/3900;
590
591 // The differential resistance of the G-APD, i.e. the dependence of the
592 // current above the breakdown voltage, is given by
593 //const double Rapd = Uov/Iapd;
594 // This allows us to estimate the current Iov at the overvoltage we want to apply
595 //const double Iov = overvoltage/Rapd;
596
597 // Estimate set point for over-voltage (voltage drop at the target point)
598 //const double Uset = Ubd + overvoltage + R*Iov*N;
599 //const double Uset = Uov<0.3 ? Ubd + overvoltage + Udrp : Ubd + overvoltage + Udrp*pow(overvoltage/Uov, 1.66);
600 const double Uset = Uov<0 ?
601 Ubd + overvoltage + Udrp*pow(overvoltage/0.34+1, 1.66) :
602 Ubd + overvoltage + Udrp*pow((overvoltage+0.34)/(Uov+0.34), 1.66);
603
604 if (fabs(overvoltage-Uov)>0.033)
605 Ndev[0]++;
606 if (fabs(overvoltage-Uov)>0.022)
607 Ndev[1]++;
608 if (fabs(overvoltage-Uov)>0.011)
609 Ndev[2]++;
610
611 // Voltage set point
612 vec[i] = Uset;
613
614 /*
615 if (fDimBias.state()==BIAS::State::kVoltageOn && GetCurrentState()==Feedback::State::kInProgress &&
616 fabs(Uov-overvoltage)>0.033)
617 cout << setprecision(4) << setw(3) << i << ": Uov=" << Uov << " Udrp=" << Udrp << " Iapd=" << Iapd*1e6 << endl;
618 */
619
620 // Calculate statistics only for channels with a valid calibration
621 //if (Uov>0)
622 {
623 const int g = hv.group();
624
625 med[g][num[g]] = Uov;
626 avg[g] += Uov;
627 num[g]++;
628
629 if (Uov<min[g])
630 min[g] = Uov;
631 if (Uov>max[g])
632 max[g] = Uov;
633
634 const double iapd = Iapd*1e6; // A --> uA
635
636 data.I[i] = iapd;
637 data.Iavg += iapd;
638 data.Irms += iapd*iapd;
639
640 data.Uov[i] = Uov;
641
642 med[2][num[2]++] = iapd;
643 }
644 }
645
646 // ------------------------------- Update voltages ------------------------------------
647
648 if (GetCurrentState()!=Feedback::State::kCalibrated) // WaitingForData, InProgress
649 {
650 if (fDimBias.state()!=BIAS::State::kRamping)
651 {
652 DimClient::sendCommandNB("BIAS_CONTROL/SET_ALL_CHANNELS_VOLTAGE",
653 vec.data(), BIAS::kNumChannels*sizeof(float));
654
655 ostringstream msg;
656 msg << setprecision(4) << "Sending new absolute offset: dU(" << fTemp << "degC)=" << fTempOffsetAvg << "V+-" << fTempOffsetRms << ", Unom=" << overvoltage << "V, Uov=" << (num[0]+num[1]>0?(avg[0]+avg[1])/(num[0]+num[1]):0) << " [N=" << Ndev[0] << "/" << Ndev[1] << "/" << Ndev[2] << "]";
657 Info(msg);
658 }
659 }
660 else
661 {
662 if (fDimBias.state()==BIAS::State::kVoltageOn)
663 {
664 ostringstream msg;
665 msg << setprecision(4) << "Current status: dU(" << fTemp << "degC)=" << fTempOffsetAvg << "V+-" << fTempOffsetRms << ", Unom=" << overvoltage << "V, Uov=" << (num[0]+num[1]>0?(avg[0]+avg[1])/(num[0]+num[1]):0) << " [N=" << Ndev[0] << "/" << Ndev[1] << "/" << Ndev[2] << "]";
666 Info(msg);
667 }
668
669 }
670
671 if (GetCurrentState()==Feedback::State::kInProgress &&
672 fDimBias.state()==BIAS::State::kRamping)
673 return GetCurrentState();
674
675 // --------------------------------- Console out --------------------------------------
676
677 if (num[0]>0 && num[1]>0 && fIsVerbose && !fDimBias.state()==BIAS::State::kRamping)
678 {
679 sort(med[0].begin(), med[0].begin()+num[0]);
680 sort(med[1].begin(), med[1].begin()+num[1]);
681
682 ostringstream msg;
683 msg << " Avg0=" << setw(7) << avg[0]/num[0] << " | Avg1=" << setw(7) << avg[1]/num[1];
684 Debug(msg);
685
686 msg.str("");
687 msg << " Med0=" << setw(7) << med[0][num[0]/2] << " | Med1=" << setw(7) << med[1][num[1]/2];
688 Debug(msg);
689
690 msg.str("");
691 msg << " Min0=" << setw(7) << min[0] << " | Min1=" << setw(7) << min[1];
692 Debug(msg);
693
694 msg.str("");
695 msg << " Max0=" << setw(7) << max[0] << " | Max1=" << setw(7) << max[1];
696 Debug(msg);
697 }
698
699 // ---------------------------- Calibrated Currents -----------------------------------
700
701 if (num[2]>0)
702 {
703 data.Iavg /= num[2];
704 data.Irms /= num[2];
705 data.Irms -= data.Iavg*data.Iavg;
706
707 data.N = num[2];
708 data.Irms = data.Irms<0 ? 0: sqrt(data.Irms);
709
710 sort(med[2].data(), med[2].data()+num[2]);
711
712 data.Imed = num[2]%2 ? med[2][num[2]/2] : (med[2][num[2]/2-1]+med[2][num[2]/2])/2;
713
714 for (int i=0; i<num[2]; i++)
715 med[2][i] = fabs(med[2][i]-data.Imed);
716
717 sort(med[2].data(), med[2].data()+num[2]);
718
719 data.Idev = med[2][uint32_t(0.682689477208650697*num[2])];
720
721 data.Tdiff = evt.GetTime().UnixTime()-fTimeCalib.UnixTime();
722
723 // FIXME:
724 // + Current overvoltage
725 // + Temp offset
726 // + User offset
727 // + Command overvoltage
728 fDimCurrents.setQuality(GetCurrentState());
729 fDimCurrents.setData(&data, sizeof(dim_data));
730 fDimCurrents.Update(evt.GetTime());
731 }
732
733 return GetCurrentState()==Feedback::State::kCalibrated ? Feedback::State::kCalibrated : Feedback::State::kInProgress;
734 }
735
736 // ======================================================================
737
738 int Print() const
739 {
740 Out() << fDim << endl;
741 Out() << fDimFSC << endl;
742 Out() << fDimBias << endl;
743
744 return GetCurrentState();
745 }
746
747 int PrintCalibration()
748 {
749 /*
750 if (fCalibration.size()==0)
751 {
752 Out() << "No calibration performed so far." << endl;
753 return GetCurrentState();
754 }
755
756 const float *avg = fCalibration.data();
757 const float *rms = fCalibration.data()+BIAS::kNumChannels;
758 const float *res = fCalibration.data()+BIAS::kNumChannels*2;
759
760 Out() << "Average current at " << fCalibrationOffset << "V below G-APD operation voltage:\n";
761
762 for (int k=0; k<13; k++)
763 for (int j=0; j<8; j++)
764 {
765 Out() << setw(2) << k << "|" << setw(2) << j*4 << "|";
766 for (int i=0; i<4; i++)
767 Out() << Tools::Form(" %6.1f+-%4.1f", avg[k*32+j*4+i], rms[k*32+j*4+i]);
768 Out() << '\n';
769 }
770 Out() << '\n';
771
772 Out() << "Measured calibration resistor:\n";
773 for (int k=0; k<13; k++)
774 for (int j=0; j<4; j++)
775 {
776 Out() << setw(2) << k << "|" << setw(2) << j*8 << "|";
777 for (int i=0; i<8; i++)
778 Out() << Tools::Form(" %5.0f", res[k*32+j*8+i]);
779 Out() << '\n';
780 }
781
782 Out() << flush;
783 */
784 return GetCurrentState();
785 }
786
787 int SetVerbosity(const EventImp &evt)
788 {
789 if (!CheckEventSize(evt.GetSize(), "SetVerbosity", 1))
790 return kSM_FatalError;
791
792 fIsVerbose = evt.GetBool();
793
794 return GetCurrentState();
795 }
796
797 int SetCurrentRequestInterval(const EventImp &evt)
798 {
799 if (!CheckEventSize(evt.GetSize(), "SetCurrentRequestInterval", 2))
800 return kSM_FatalError;
801
802 fCurrentRequestInterval = evt.GetUShort();
803
804 Info("New current request interval: "+to_string(fCurrentRequestInterval)+"ms");
805
806 return GetCurrentState();
807 }
808
809 int Calibrate()
810 {
811 if (fDimBias.state()!=BIAS::State::kVoltageOff)
812 {
813 Warn("Calibration can only be started when biasctrl is in state VoltageOff.");
814 return GetCurrentState();
815 }
816
817 Message("Starting calibration (ignore="+to_string(fNumCalibIgnore)+", N="+to_string(fNumCalibRequests)+")");
818
819 fCursorCur = -fNumCalibIgnore;
820 fCurrentsAvg.assign(BIAS::kNumChannels, 0);
821 fCurrentsRms.assign(BIAS::kNumChannels, 0);
822
823 fBiasDac.assign(BIAS::kNumChannels, 0);
824
825 fCalibStep = 3;
826 fTimeCalib = Time();
827
828 Dim::SendCommandNB("BIAS_CONTROL/SET_GLOBAL_DAC", uint32_t(256+512*fCalibStep));
829
830 return Feedback::State::kCalibrating;
831 }
832
833 int Start(const EventImp &evt)
834 {
835 if (!CheckEventSize(evt.GetSize(), "Start", 4))
836 return kSM_FatalError;
837
838 if (fDimBias.state()==BIAS::State::kRamping)
839 {
840 Warn("Feedback can not be started when biasctrl is in state Ramping.");
841 return GetCurrentState();
842 }
843
844 fUserOffset = evt.GetFloat();
845
846 fCursorCur = 0;
847
848 fCurrentsAvg.assign(BIAS::kNumChannels, 0);
849 fCurrentsRms.assign(BIAS::kNumChannels, 0);
850
851 ostringstream out;
852 out << "Starting feedback with an offset of " << fUserOffset << "V";
853 Message(out);
854
855 return Feedback::State::kWaitingForData;
856 }
857
858 int StopFeedback()
859 {
860 if (GetCurrentState()==Feedback::State::kCalibrating)
861 return Feedback::State::kConnected;
862
863 if (GetCurrentState()>Feedback::State::kCalibrated)
864 return Feedback::State::kCalibrated;
865
866 return GetCurrentState();
867 }
868
869 bool LoadOffsets(const string &file)
870 {
871 vector<double> data(416);
872
873 ifstream fin(file);
874
875 int cnt = 0;
876 while (fin && cnt<320)
877 fin >> data[cnt++];
878
879 if (cnt!=320)
880 {
881 Error("Reading offsets from "+file+" failed [N="+to_string(cnt-1)+"]");
882 return false;
883 }
884
885 fVoltOffset = data;
886
887 fDimOffsets.Update(fVoltOffset);
888
889 Info("New voltage offsets loaded from "+file);
890 return true;
891
892 }
893
894 int LoadOffset(const EventImp &evt)
895 {
896 LoadOffsets(evt.GetText());
897 return GetCurrentState();
898 }
899
900 int ResetOffset()
901 {
902 fVoltOffset.assign(416, 0);
903
904 fDimOffsets.Update(fVoltOffset);
905
906 Info("Voltage offsets resetted.");
907 return GetCurrentState();
908 }
909
910
911 int Execute()
912 {
913 if (!fDim.online())
914 return Feedback::State::kDimNetworkNA;
915
916 const bool bias = fDimBias.state() >= BIAS::State::kConnecting;
917 const bool fsc = fDimFSC.state() >= FSC::State::kConnected;
918
919 // All subsystems are not connected
920 if (!bias && !fsc)
921 return Feedback::State::kDisconnected;
922
923 // Not all subsystems are yet connected
924 if (!bias || !fsc)
925 return Feedback::State::kConnecting;
926
927 if (GetCurrentState()<Feedback::State::kCalibrating)
928 return Feedback::State::kConnected;
929
930 if (GetCurrentState()==Feedback::State::kConnected)
931 return GetCurrentState();
932 if (GetCurrentState()==Feedback::State::kCalibrating)
933 return GetCurrentState();
934
935 // kCalibrated, kWaitingForData, kInProgress
936
937 if (fDimBias.state()==BIAS::State::kVoltageOn || (fDimBias.state()==BIAS::State::kVoltageOff && GetCurrentState()==Feedback::State::kWaitingForData))
938 {
939 static Time past;
940 if (fCurrentRequestInterval>0 && Time()-past>boost::posix_time::milliseconds(fCurrentRequestInterval))
941 {
942 Dim::SendCommandNB("BIAS_CONTROL/REQUEST_STATUS");
943 past = Time();
944 }
945 }
946
947 return GetCurrentState();
948 }
949
950public:
951 StateMachineFeedback(ostream &out=cout) : StateMachineDim(out, "FEEDBACK"),
952 fIsVerbose(false),
953 //---
954 fDimFSC("FSC_CONTROL"),
955 fDimBias("BIAS_CONTROL"),
956 //---
957 fDimCalibration("FEEDBACK/CALIBRATION", "F:416;F:416;F:416;F:416",
958 "Current offsets"
959 "|Avg[uA]:Average offset at dac=256+5*512"
960 "|Rms[uA]:Rms of Avg"
961 "|R[Ohm]:Measured calibration resistor"
962 "|U[V]:Corresponding voltage reported by biasctrl"),
963 fDimCalibration2("FEEDBACK/CALIBRATION_STEPS", "I:1;F:416;F:416;F:416",
964 "Calibration of the R8 resistor"
965 "|DAC[dac]:DAC setting"
966 "|U[V]:Corresponding voltages reported by biasctrl"
967 "|Iavg[uA]:Averaged measured current"
968 "|Irms[uA]:Rms measured current"),
969 fDimCalibrationR8("FEEDBACK/CALIBRATION_R8", "F:416;F:416",
970 "Calibration of R8"
971 "|DeltaI[uA]:Average offset"
972 "|R8[uA]:Measured effective resistor R8"),
973 fDimCurrents("FEEDBACK/CALIBRATED_CURRENTS", "F:416;F:1;F:1;F:1;F:1;I:1;F:1;F:416;F:1;F:1",
974 "Calibrated currents"
975 "|I[uA]:Calibrated currents per pixel"
976 "|I_avg[uA]:Average calibrated current (N channels)"
977 "|I_rms[uA]:Rms of calibrated current (N channels)"
978 "|I_med[uA]:Median calibrated current (N channels)"
979 "|I_dev[uA]:Deviation of calibrated current (N channels)"
980 "|N[uint16]:Number of valid values"
981 "|T_diff[s]:Time difference to calibration"
982 "|U_ov[V]:Calculated overvoltage"
983 "|U_nom[V]:Nominal overvoltage"
984 "|dU_temp[V]:Correction calculated from temperature"
985 ),
986 fDimOffsets("FEEDBACK/OFFSETS", "F:416",
987 "Offsets operation voltages"
988 "|U[V]:Offset per bias channels"),
989 fVoltOffset(416),
990 fCurrentRequestInterval(0),
991 fNumCalibIgnore(30),
992 fNumCalibRequests(300)
993 {
994 fDim.Subscribe(*this);
995 fDimFSC.Subscribe(*this);
996 fDimBias.Subscribe(*this);
997
998 fDimBias.SetCallback(bind(&StateMachineFeedback::HandleBiasStateChange, this));
999
1000 Subscribe("BIAS_CONTROL/CURRENT")
1001 (bind(&StateMachineFeedback::HandleBiasCurrent, this, placeholders::_1));
1002 Subscribe("BIAS_CONTROL/VOLTAGE")
1003 (bind(&StateMachineFeedback::HandleBiasVoltage, this, placeholders::_1));
1004 Subscribe("BIAS_CONTROL/DAC")
1005 (bind(&StateMachineFeedback::HandleBiasDac, this, placeholders::_1));
1006 Subscribe("BIAS_CONTROL/NOMINAL")
1007 (bind(&StateMachineFeedback::HandleBiasNom, this, placeholders::_1));
1008 Subscribe("FSC_CONTROL/BIAS_TEMP")
1009 (bind(&StateMachineFeedback::HandleCameraTemp, this, placeholders::_1));
1010
1011 // State names
1012 AddStateName(Feedback::State::kDimNetworkNA, "DimNetworkNotAvailable",
1013 "The Dim DNS is not reachable.");
1014
1015 AddStateName(Feedback::State::kDisconnected, "Disconnected",
1016 "The Dim DNS is reachable, but the required subsystems are not available.");
1017 AddStateName(Feedback::State::kConnecting, "Connecting",
1018 "Either biasctrl or fscctrl not connected.");
1019 AddStateName(Feedback::State::kConnected, "Connected",
1020 "biasctrl and fscctrl are available and connected with their hardware.");
1021
1022 AddStateName(Feedback::State::kCalibrating, "Calibrating",
1023 "Bias crate calibrating in progress.");
1024 AddStateName(Feedback::State::kCalibrated, "Calibrated",
1025 "Bias crate calibrated.");
1026
1027 AddStateName(Feedback::State::kWaitingForData, "WaitingForData",
1028 "Current control started, waiting for valid temperature and current data.");
1029 AddStateName(Feedback::State::kInProgress, "InProgress",
1030 "Current control in progress.");
1031
1032
1033 /*
1034 AddEvent("SET_CURRENT_REQUEST_INTERVAL")
1035 (bind(&StateMachineFeedback::SetCurrentRequestInterval, this, placeholders::_1))
1036 ("|interval[ms]:Interval between two current requests in modes which need that.");
1037 */
1038
1039 AddEvent("CALIBRATE", Feedback::State::kConnected, Feedback::State::kCalibrated)
1040 (bind(&StateMachineFeedback::Calibrate, this))
1041 ("");
1042
1043 AddEvent("START", "F:1", Feedback::State::kCalibrated)
1044 (bind(&StateMachineFeedback::Start, this, placeholders::_1))
1045 ("Start the current/temperature control loop"
1046 "|Uov[V]:Overvoltage to be applied (standard value is 1.1V)");
1047
1048 AddEvent("STOP")
1049 (bind(&StateMachineFeedback::StopFeedback, this))
1050 ("Stop any control loop");
1051
1052 AddEvent("LOAD_OFFSETS", "C", Feedback::State::kConnected, Feedback::State::kCalibrated)
1053 (bind(&StateMachineFeedback::LoadOffset, this, placeholders::_1))
1054 ("");
1055 AddEvent("RESET_OFFSETS", Feedback::State::kConnected, Feedback::State::kCalibrated)
1056 (bind(&StateMachineFeedback::ResetOffset, this))
1057 ("");
1058
1059
1060 AddEvent("PRINT")
1061 (bind(&StateMachineFeedback::Print, this))
1062 ("");
1063 AddEvent("PRINT_CALIBRATION")
1064 (bind(&StateMachineFeedback::PrintCalibration, this))
1065 ("");
1066
1067 // Verbosity commands
1068 AddEvent("SET_VERBOSE", "B:1")
1069 (bind(&StateMachineFeedback::SetVerbosity, this, placeholders::_1))
1070 ("set verbosity state"
1071 "|verbosity[bool]:disable or enable verbosity when calculating overvoltage");
1072 }
1073
1074 int EvalOptions(Configuration &conf)
1075 {
1076 fIsVerbose = !conf.Get<bool>("quiet");
1077
1078 if (!fMap.Read(conf.Get<string>("pixel-map-file")))
1079 {
1080 Error("Reading mapping table from "+conf.Get<string>("pixel-map-file")+" failed.");
1081 return 1;
1082 }
1083
1084 fCurrentRequestInterval = conf.Get<uint16_t>("current-request-interval");
1085 fNumCalibIgnore = conf.Get<uint16_t>("num-calib-ignore");
1086 fNumCalibRequests = conf.Get<uint16_t>("num-calib-average");
1087 fTempCoefficient = conf.Get<double>("temp-coefficient");
1088
1089 if (conf.Has("offset-file"))
1090 LoadOffsets(conf.Get<string>("offset-file"));
1091
1092 return -1;
1093 }
1094};
1095
1096// ------------------------------------------------------------------------
1097
1098#include "Main.h"
1099
1100template<class T>
1101int RunShell(Configuration &conf)
1102{
1103 return Main::execute<T, StateMachineFeedback>(conf);
1104}
1105
1106void SetupConfiguration(Configuration &conf)
1107{
1108 po::options_description control("Feedback options");
1109 control.add_options()
1110 ("quiet,q", po_bool(true), "Disable printing more information on average overvoltagecontents of all received messages (except dynamic data) in clear text.")
1111 ("pixel-map-file", var<string>()->required(), "Pixel mapping file. Used here to get the default reference voltage.")
1112 ("current-request-interval", var<uint16_t>(1000), "Interval between two current requests.")
1113 ("num-calib-ignore", var<uint16_t>(30), "Number of current requests to be ignored before averaging")
1114 ("num-calib-average", var<uint16_t>(300), "Number of current requests to be averaged")
1115 ("temp-coefficient", var<double>()->required(), "Temp. coefficient [V/K]")
1116 ("offset-file", var<string>(), "File with operation voltage offsets")
1117 ;
1118
1119 conf.AddOptions(control);
1120}
1121
1122/*
1123 Extract usage clause(s) [if any] for SYNOPSIS.
1124 Translators: "Usage" and "or" here are patterns (regular expressions) which
1125 are used to match the usage synopsis in program output. An example from cp
1126 (GNU coreutils) which contains both strings:
1127 Usage: cp [OPTION]... [-T] SOURCE DEST
1128 or: cp [OPTION]... SOURCE... DIRECTORY
1129 or: cp [OPTION]... -t DIRECTORY SOURCE...
1130 */
1131void PrintUsage()
1132{
1133 cout <<
1134 "The feedback control the BIAS voltages based on the calibration signal.\n"
1135 "\n"
1136 "The default is that the program is started without user intercation. "
1137 "All actions are supposed to arrive as DimCommands. Using the -c "
1138 "option, a local shell can be initialized. With h or help a short "
1139 "help message about the usuage can be brought to the screen.\n"
1140 "\n"
1141 "Usage: feedback [-c type] [OPTIONS]\n"
1142 " or: feedback [OPTIONS]\n";
1143 cout << endl;
1144}
1145
1146void PrintHelp()
1147{
1148 Main::PrintHelp<StateMachineFeedback>();
1149
1150 /* Additional help text which is printed after the configuration
1151 options goes here */
1152
1153 /*
1154 cout << "bla bla bla" << endl << endl;
1155 cout << endl;
1156 cout << "Environment:" << endl;
1157 cout << "environment" << endl;
1158 cout << endl;
1159 cout << "Examples:" << endl;
1160 cout << "test exam" << endl;
1161 cout << endl;
1162 cout << "Files:" << endl;
1163 cout << "files" << endl;
1164 cout << endl;
1165 */
1166}
1167
1168int main(int argc, const char* argv[])
1169{
1170 Configuration conf(argv[0]);
1171 conf.SetPrintUsage(PrintUsage);
1172 Main::SetupConfiguration(conf);
1173 SetupConfiguration(conf);
1174
1175 if (!conf.DoParse(argc, argv, PrintHelp))
1176 return 127;
1177
1178 //try
1179 {
1180 // No console access at all
1181 if (!conf.Has("console"))
1182 {
1183// if (conf.Get<bool>("no-dim"))
1184// return RunShell<LocalStream, StateMachine, ConnectionFSC>(conf);
1185// else
1186 return RunShell<LocalStream>(conf);
1187 }
1188 // Cosole access w/ and w/o Dim
1189/* if (conf.Get<bool>("no-dim"))
1190 {
1191 if (conf.Get<int>("console")==0)
1192 return RunShell<LocalShell, StateMachine, ConnectionFSC>(conf);
1193 else
1194 return RunShell<LocalConsole, StateMachine, ConnectionFSC>(conf);
1195 }
1196 else
1197*/ {
1198 if (conf.Get<int>("console")==0)
1199 return RunShell<LocalShell>(conf);
1200 else
1201 return RunShell<LocalConsole>(conf);
1202 }
1203 }
1204 /*catch (std::exception& e)
1205 {
1206 cerr << "Exception: " << e.what() << endl;
1207 return -1;
1208 }*/
1209
1210 return 0;
1211}
Note: See TracBrowser for help on using the repository browser.