Matlab: polyfit.m void GammaAnalyALG::calEnergyUpdate1(uvec idx, mat ELib) { // get pat centroids vec C = pat.Peaks.col(i_Centroid); C = C(idx); // fit parameter //[s, msg, p, perr] = calFitPara('Energy', C, ELib, []); mat p, perr, ep; bool s = calFitPara(p, perr, Cal_Energy, C, ELib); if (s) // success, set parameter, data points and flag { //[s, msg, name] = setCalibrationData(sn, 'Energy', 'CalUpdate1', C, ELib, []); QString name, pname; s = setCalibrationData(name, Cal_Energy, "CalUpdate1", C, ELib); if(s) { //[s, pname] = setCalibrationPara(sn, 'Energy', 'CalUpdate1', p, [], name); s = setCalibrationPara(pname, Cal_Energy, "CalUpdate1", p, ep, name); //lPrintSetFail(sn, s, msg); if(!s) { qDebug() << "Failed saving energy calibration parameter in calEnergyUpdate1"; } } else { qDebug() << "Failed saving energy calibration data points"; //[s, pname] = setCalibrationPara(sn, 'Energy', 'CalUpdate1', p, [], ''); s = setCalibrationPara(pname, Cal_Energy, "CalUpdate1", p, ep, ""); //lPrintSetFail(sn, s, msg); if(!s) { qDebug() << "Failed saving energy calibration parameter in calEnergyUpdate1"; } } //if (s) { s = setCurrentCalibration(Cal_Energy, pname); } } else { qDebug() << QString("First energy calibration update failed!"); } } void GammaAnalyALG::calResUpdate(uvec idx) { // Peaks required for fitting int NPeaksMin = 4; uword NPeaks = idx.size(); //[NA, C, FC] = dmps(sn, 'NetArea', 'Centroid', 'FWHM_Ch'); QStringList Opts; Opts << "NetArea" << "Centroid" << "FWHM_Ch"; field t_f = dmps(Opts); vec NA = t_f(0)(idx); vec C = t_f(1)(idx); vec FC = t_f(2)(idx); // fit peaks uvec Bool = zeros(NPeaks); uvec tmp; vec ret1, ret2, ret3; double X2; for(uword i=0; i NPeaksMin) { // fit new resolution //[s, msg, p, perr] = calFitPara('Resolution', E(idx), FE(idx), []); s = calFitPara(p, perr, Cal_Resolution, E(idx), FE(idx)); } else { qDebug() << "Not enough peaks for fitting"; s = false; } // write results to global if (s) { // success, set parameters QString pname; if (stat_sv) { //[s, pname] = setCalibrationPara(sn, 'Resolution', 'CalUpdate', p, perr, name); s = setCalibrationPara(pname, Cal_Resolution, "CalUpdate", p, perr, name); //lPrintSetFail(sn, s, msg) } else { //[s, pname] = setCalibrationPara(sn, 'Resolution', 'CalUpdate', p, perr, ''); s = setCalibrationPara(pname, Cal_Resolution, "CalUpdate", p, perr, ""); //lPrintSetFail(sn, s, msg) } if (s) { //setCurrentCalibration(sn, 'Resolution', pname); setCurrentCalibration(Cal_Resolution, pname); } } else { qDebug() << "Resolution Calibration Update Failed"; } } void GammaAnalyALG::calEnergyUpdate2(uvec idx, mat ELib) { uword NPeaks = idx.size(); vec C = zeros(NPeaks); // fit peaks vec ret1, ret2, ret3; uvec tmp; double X2; bool s; for(uword i=0; i 99) return false; name = QString("%1 %2").arg(src).arg(num); // Check if new data is valid, reshape if necessary mat data; bRet = lCalDataCheck(data, x, y, err); if (!bRet) { qDebug() << QString("Invalid calibration data from source %1").arg(src); return bRet; } // pack and write list back to global /*caldatalist(pos, 1) = {name}; caldatalist(pos, 2) = {data}; s = setCalDataGlob(sn, cal, caldatalist);*/ QMap& subMap = m_CalData[cal]; subMap[name] = data; return bRet; } bool GammaAnalyALG::setCalibrationPara(QString &pname, CalibrationType cal, QString src, mat p, mat perr, QString dnm) { bool bRet = true; if(p.is_empty()) return false; /*[s, msg, calentries] = getCalParaGlob(sn, cal); if ~s logWrite('%s', msg); end names = calentries(:, 1); [s, msg, name, pos] = lGetNewName(src, names); if ~s logWrite('invalid %s calibration source for spectrum nr %d: %s', ... cal, sn, msg); return end*/ QMap2Para::iterator it = m_CalPara.find(cal); QStringList names; int num = 1; if(it != m_CalPara.end()) { names = it.value().keys(); num = names.filter(src).size() + 1; } if(src.contains(' ') || num > 99) return false; pname = QString("%1 %2").arg(src).arg(num); // handle interpolation: internal p also contains x,y data if(p.size() == 1 && p(0u) == 1) //if isequal(p, 1) { bRet = lMakeInterp(p, perr, cal, m_CalData[cal][dnm]); if (!bRet) { return bRet; } } //[s, msg] = lCalParaCheck(cal, p); bRet = Independ::calParaCheck(p(0u), p.cols(1u, p.size()-1u)); if (!bRet) return bRet; //[s, msg, perr] = lCalCheckPerr(p, perr); if (perr.is_empty()) perr = gNaN * ones(p.size()-1u, 1u); else { if (perr.size() != p.size() - 1u) { qDebug() << "error estimate of wrong size"; return false; } else if (any(vectorise(perr) < 0)) { qDebug() << "negative error estimate"; return false; } } /*if isempty(data) dnm = ''; elseif isstr(data) [s, msg, x, y, err, dnm] = getCalibrationData(sn, cal, data); else [s, msg, dnm] = setCalibrationData(sn, cal, src, data); if ~s logWrite(... 'failed storing data set, calib [%s], source [%s], sn [%d]: %s', ... cal, src, sn, msg); return end end [s, msg] = lCalCheckData(sn, cal, dnm); if ~s logWrite( ... 'invalid data name, calib [%s], source [%s], sn [%d]: %s', ... cal, src, sn, msg); return end*/ /*calentries(pos, 1) = {name}; calentries(pos, 2) = {p}; calentries(pos, 3) = {perr}; calentries(pos, 4) = {dnm}; s = setCalParaGlob(sn, cal, calentries);*/ CalibPara cp; cp.name = pname; cp.dnm = dnm; cp.p = p; cp.perr = perr; QMap& subMap = m_CalPara[cal]; subMap[pname] = cp; return bRet; } bool GammaAnalyALG::setCurrentCalibration(CalibrationType cal, QString name) { // check only one spectrum is affected bool s = true; // get index in list of calibration names vec sqi, uqi; field f; f << vec("1"); int n; switch (cal) { case Cal_Energy: n = 1; break; case Cal_Resolution: patSetPeaks(uvec(0u), QStringList("ResChanged"), f, true); n = 2; break; case Cal_Efficiency: n = 3; break; case Cal_Tot_efficiency: n = 4; break; case Cal_Tail: patSetPeaks(uvec(0u), QStringList("TailChanged"), f, true); n = 5; break; case Cal_Tail_alpha: n = 6; break; case Cal_Tail_right: n = 7; break; case Cal_Tail_right_alpha: n = 8; break; case Cal_Step_ratio: patBaseVar(sqi, uqi); n = 9; break; default: s = false; qDebug() << "invalid calibration string in setCurrentCalibration"; } // get old PAT /*[s, msg, patEntry] = getPATByName(sn, 'TrueCurrent'); patName = patEntry{1}; calData = patEntry{5}; nat = patEntry{7}; if ~s return end // missing: check that name is valid calibration string calData{n} = name; // raise update-flag in NAT nat{6} = 1; // write change to global patEntry([5,7]) = {calData, nat}; setPATEntry(sn, patEntry); // change in temporary as well, if there is a temporary pat [s, msg, tpatEntry] = getPATByName(sn, 'Temporary'); tpatName = tpatEntry{1}; tcalData = tpatEntry{5}; tnat = tpatEntry{7}; if ~s | isempty(tpatName) return end tcalData{n} = name; tnat{6} = 1; tpatEntry([5,7]) = {tcalData, tnat}; setPATEntry(sn, tpatEntry);*/ pat.CalibName[n-1] = name; }