#include #include "Dim.h" #include "Event.h" #include "Shell.h" #include "StateMachineDim.h" #include "Connection.h" #include "Configuration.h" #include "Console.h" #include "Converter.h" #include "DimServiceInfoList.h" #include "PixelMap.h" #include "tools.h" #include "LocalControl.h" #include "HeadersFAD.h" #include "HeadersBIAS.h" namespace ba = boost::asio; namespace bs = boost::system; namespace dummy = ba::placeholders; using namespace std; // ------------------------------------------------------------------------ #include "DimDescriptionService.h" // ------------------------------------------------------------------------ class StateMachineFeedback : public StateMachineDim, public DimInfoHandler { /* int Wrap(boost::function f) { f(); return T::GetCurrentState(); } boost::function Wrapper(boost::function func) { return bind(&StateMachineMCP::Wrap, this, func); }*/ private: enum states_t { kStateDimNetworkNA = 1, kStateDisconnected, kStateConnecting, kStateConnected, kStateTempCtrlIdle, kStateFeedbackCtrlIdle, kStateTempCtrlRunning, kStateFeedbackCtrlRunning, }; enum control_t { kIdle, kTemp, kFeedback }; control_t fControlType; PixelMap fMap; DimServiceInfoList fNetwork; pair fStatusDim; pair fStatusFAD; pair fStatusFSC; pair fStatusBias; DimStampedInfo fDim; DimStampedInfo fFAD; DimStampedInfo fFSC; DimStampedInfo fBias; DimStampedInfo fBiasData; DimStampedInfo fCameraTemp; DimDescribedService fDimReference; DimDescribedService fDimDeviation; vector> fData; uint64_t fCursor; Time fBiasLast; Time fStartTime; valarray fPV[3]; // Process variable (intgerated/averaged amplitudes) valarray fSP; // Set point (target amplitudes) double fKp; // Proportional constant double fKi; // Integral constant double fKd; // Derivative constant double fT; // Time constant (cycle time) double fGain; // Gain (conversion from a DRS voltage deviation into a BIAS voltage change at G-APD reference voltage) double fT21; double fBiasOffset; bool fOutputEnabled; pair GetNewState(DimStampedInfo &info) const { const bool disconnected = info.getSize()==0; // Make sure getTimestamp is called _before_ getTimestampMillisecs const int tsec = info.getTimestamp(); const int tms = info.getTimestampMillisecs(); return make_pair(Time(tsec, tms*1000), disconnected ? -2 : info.getQuality()); } void ResetData() { fData.clear(); fData.resize(500); fCursor = 0; fStartTime = Time(); fSP.resize(416); fSP = valarray(0., 416); vector vec(2*BIAS::kNumChannels); fDimDeviation.Update(vec); fPV[0].resize(0); fPV[1].resize(0); fPV[2].resize(0); if (fKp==0 && fKi==0 && fKd==0) Warn("Control loop parameters are all set to zero."); } void infoHandler() { DimInfo *curr = getInfo(); // get current DimInfo address if (!curr) return; if (curr==&fBias) { fStatusBias = GetNewState(fBias); return; } if (curr==&fFAD) { fStatusFAD = GetNewState(fFAD); return; } if (curr==&fFSC) { fStatusFSC = GetNewState(fFSC); return; } if (curr==&fDim) { fStatusDim = GetNewState(fDim); fStatusDim.second = curr->getSize()==4 ? curr->getInt() : 0; return; } if (curr==&fCameraTemp) { if (curr->getSize()==0) return; if (curr->getSize()!=60*sizeof(float)) return; const float *ptr = (float*)curr->getData(); double avg = 0; int num = 0; for (int i=1; i<32; i++) if (ptr[i]!=0) { avg += ptr[i]; num++; } if (num==0) return; avg /= num; const float diff = (avg-25)*4./70 + fBiasOffset; vector vec(2*BIAS::kNumChannels); for (int i=0; igetSize()==0) return; if (curr->getSize()!=1440*sizeof(float)) return; // -------- Check age of last stored event -------- // Must be called in this order const int tsec = curr->getTimestamp(); const int tms = curr->getTimestampMillisecs(); const Time tm(tsec, tms*1000); if (Time()-fBiasLast>boost::posix_time::seconds(30)) { Warn("Last received event data older than 30s... resetting average calculation."); ResetData(); } fBiasLast = tm; // -------- Store new event -------- fData[fCursor%fData.size()].assign(reinterpret_cast(curr->getData()), reinterpret_cast(curr->getData())+1440); fCursor++; if (fCursor med(1440); for (int ch=0; ch<1440; ch++) { vector arr(fData.size()); for (size_t i=0; i med(1440); vector rms(1440); for (size_t i=0; i avg(BIAS::kNumChannels); vector num(BIAS::kNumChannels); for (int i=0; i<1440; i++) { const PixelMapEntry &ch = fMap.hw(i); // FIXME: Add a consistency check if the median makes sense... // FIXME: Add a consistency check to remove pixels with bright stars (median?) avg[ch.hv()] += med[i]; num[ch.hv()]++; } for (int i=0; i0 ? fT : (tm0-fStartTime).total_microseconds()/1000000.; const double T10 = fT21; fT21 = T21; fStartTime = tm0; ostringstream out; out << "New " << fData.size() << " event received: " << fCursor << " / " << setprecision(3) << T21 << "s"; Info(out); if (fPV[0].size()==0) { fPV[0].resize(avg.size()); fPV[0] = valarray(avg.data(), avg.size()); } else if (fPV[1].size()==0) { fPV[1].resize(avg.size()); fPV[1] = valarray(avg.data(), avg.size()); } else if (fPV[2].size()==0) { fPV[2].resize(avg.size()); fPV[2] = valarray(avg.data(), avg.size()); } else { fPV[0] = fPV[1]; fPV[1] = fPV[2]; fPV[2].resize(avg.size()); fPV[2] = valarray(avg.data(), avg.size()); if (T10<=0 || T21<=0) return; cout << "Calculating (" << fCursor << ":" << T21 << ")... " << endl; // fKi[j] = response[j]*gain; // Kp = 0; // Kd = 0; // => Kp = 0.01 * gain = 0.00005 // => Ki = 0.8 * gain/20s = 0.00025 // => Kd = 0.1 * gain/20s = 0.00003 //valarray correction = // - Kp * (PV[2] - PV[1]) // + dT * Ki * (SP - PV[2]) // - Kd / dT * (PV[2] - 2*PV[1] + PV[0]); // // - (Kp+Kd/dT1) * (PV[2] - PV[1]) // + dT2 * Ki * (SP - PV[2]) // + Kd / dT1 * (PV[1] - PV[0]); // // - Kp * (PV[2] - PV[1]) // + Ki * (SP - PV[2])*dT // - Kd * (PV[2] - PV[1])/dT // + Kd * (PV[1] - PV[0])/dT; // //valarray correction = // - Kp*(PV[2] - PV[1]) + Ki * T21 * (SP-PV[2]) - Kd*(PV[2]-PV[1])/T21 - Kd*(PV[0]-PV[1])/T01; const valarray correction = fGain/1000* ( - (fKp+fKd/T21)*(fPV[2] - fPV[1]) + fKi*T21*(fSP-fPV[2]) + fKd/T10*(fPV[1]-fPV[0]) ); /* integral = 0 start: integral += (fSP - fPV[2])*dt output = Kp*(fSP - fPV[2]) + Ki*integral - Kd*(fPV[2] - fPV[1])/dt wait(dt) goto start */ vector vec(2*BIAS::kNumChannels); for (int i=0; i correction(416); vector response(416); vector dev(416); double gain = 0; for (int j=0; j<416; j++) { //avg[j] /= fData.size(); //rms[j] /= sqrt((rms[j]*fData.size() - avg[j]*avg[j]))/fData.size(); dev[j] = avg[j] - target[j]; // Determine correction from response maxtrix and change voltages correction[j] = -dev[j]*response[j]*gain; if (correction[j]==0) continue; // Limit voltage steps // Limit changes to 100 mV // if (fabs(correction[j]) > 0.1) // correction[j] = 0.1*fabs(correction[j])/correction[j]; // Add voltage change command if not too noisy // if (fabs(Average[i]) > 2*Sigma[i]) // Cmd << fIDTable[i] << " " << std::showpos << Correction/Multiplicity[i] << " "; } */ return; } } bool CheckEventSize(size_t has, const char *name, size_t size) { if (has==size) return true; ostringstream msg; msg << name << " - Received event has " << has << " bytes, but expected " << size << "."; Fatal(msg); return false; } void PrintState(const pair &state, const char *server) { const State rc = fNetwork.GetState(server, state.second); Out() << state.first.GetAsStr("%H:%M:%S.%f").substr(0, 12) << " - "; Out() << kBold << server << ": "; Out() << rc.name << "[" << rc.index << "]"; Out() << kReset << " - " << kBlue << rc.comment << endl; } int Print() { Out() << fStatusDim.first.GetAsStr("%H:%M:%S.%f").substr(0, 12) << " - "; Out() << kBold << "DIM_DNS: "; if (fStatusDim.second==0) Out() << "Offline" << endl; else Out() << "V" << fStatusDim.second/100 << 'r' << fStatusDim.second%100 << endl; PrintState(fStatusFAD, "FAD_CONTROL"); PrintState(fStatusBias, "BIAS_CONTROL"); return GetCurrentState(); } int SetConstant(const EventImp &evt, int constant) { if (!CheckEventSize(evt.GetSize(), "SetConstant", 8)) return kSM_FatalError; switch (constant) { case 0: fKi = evt.GetDouble(); break; case 1: fKp = evt.GetDouble(); break; case 2: fKd = evt.GetDouble(); break; case 3: fT = evt.GetDouble(); break; case 4: fGain = evt.GetDouble(); break; default: Fatal("SetConstant got an unexpected constant id -- this is a program bug!"); return kSM_FatalError; } return GetCurrentState(); } int EnableOutput(const EventImp &evt) { if (!CheckEventSize(evt.GetSize(), "EnableOutput", 1)) return kSM_FatalError; fOutputEnabled = evt.GetBool(); return GetCurrentState(); } int StartFeedback() { ResetData(); fControlType = kFeedback; return GetCurrentState(); } int StartTempCtrl(const EventImp &evt) { if (!CheckEventSize(evt.GetSize(), "StartTempCtrl", 4)) return kSM_FatalError; ostringstream out; out << "Starting temperature feedback with an offset of " << fBiasOffset << "V"; Message(out); fBiasOffset = evt.GetFloat(); fControlType = kTemp; return GetCurrentState(); } int StopFeedback() { fControlType = kIdle; return GetCurrentState(); } int StoreReference() { if (!fPV[0].size() && !fPV[1].size() && !fPV[2].size()) { Warn("No values in memory. Take enough events first!"); return GetCurrentState(); } // FIXME: Check age if (!fPV[1].size() && !fPV[2].size()) fSP = fPV[0]; if (!fPV[2].size()) fSP = fPV[1]; else fSP = fPV[2]; vector vec(BIAS::kNumChannels); for (int i=0; i vec(BIAS::kNumChannels); for (int i=0; i("pixel-map-file"))) { Error("Reading mapping table from "+conf.Get("pixel-map-file")+" failed."); return 1; } // -110 / -110 (-23 DAC / -0.51V) // Reference voltage: -238 / -203 // -360 / -343 ( 23 DAC / 0.51V) // 0.005 A/V // 220 Amplitude / 1V // Gain = 1V / 200 = 0.005 fGain = 5; // (BIAS)V / (DRS)V ( 1V / 0.22V ) fKp = 0; fKd = 0; fKi = 0.12; fT = 1; ostringstream msg; msg << "Control loop parameters: "; msg << "Kp=" << fKp << ", Kd=" << fKd << ", Ki=" << fKi << ", "; if (fT>0) msg << fT; else msg << ""; msg << ", Gain(BIAS/DRS)=" << fGain << "V/V"; Message(msg); return -1; } }; // ------------------------------------------------------------------------ #include "Main.h" template int RunShell(Configuration &conf) { return Main::execute(conf); } void SetupConfiguration(Configuration &conf) { po::options_description control("Feedback options"); control.add_options() ("pixel-map-file", var("FACTmapV5a.txt"), "Pixel mapping file. Used here to get the default reference voltage.") ; conf.AddOptions(control); } /* Extract usage clause(s) [if any] for SYNOPSIS. Translators: "Usage" and "or" here are patterns (regular expressions) which are used to match the usage synopsis in program output. An example from cp (GNU coreutils) which contains both strings: Usage: cp [OPTION]... [-T] SOURCE DEST or: cp [OPTION]... SOURCE... DIRECTORY or: cp [OPTION]... -t DIRECTORY SOURCE... */ void PrintUsage() { cout << "The feedback control the BIAS voltages based on the calibration signal.\n" "\n" "The default is that the program is started without user intercation. " "All actions are supposed to arrive as DimCommands. Using the -c " "option, a local shell can be initialized. With h or help a short " "help message about the usuage can be brought to the screen.\n" "\n" "Usage: feedback [-c type] [OPTIONS]\n" " or: feedback [OPTIONS]\n"; cout << endl; } void PrintHelp() { Main::PrintHelp(); /* Additional help text which is printed after the configuration options goes here */ /* cout << "bla bla bla" << endl << endl; cout << endl; cout << "Environment:" << endl; cout << "environment" << endl; cout << endl; cout << "Examples:" << endl; cout << "test exam" << endl; cout << endl; cout << "Files:" << endl; cout << "files" << endl; cout << endl; */ } int main(int argc, const char* argv[]) { Configuration conf(argv[0]); conf.SetPrintUsage(PrintUsage); Main::SetupConfiguration(conf); SetupConfiguration(conf); if (!conf.DoParse(argc, argv, PrintHelp)) return -1; //try { // No console access at all if (!conf.Has("console")) { // if (conf.Get("no-dim")) // return RunShell(conf); // else return RunShell(conf); } // Cosole access w/ and w/o Dim /* if (conf.Get("no-dim")) { if (conf.Get("console")==0) return RunShell(conf); else return RunShell(conf); } else */ { if (conf.Get("console")==0) return RunShell(conf); else return RunShell(conf); } } /*catch (std::exception& e) { cerr << "Exception: " << e.what() << endl; return -1; }*/ return 0; }