27#include <boost/container/flat_set.hpp>
28#include <boost/property_tree/ptree.hpp>
29#include <boost/property_tree/xml_parser.hpp>
31#include <boost/filesystem.hpp>
32#include <boost/math/special_functions/round.hpp>
48Logger logger(
"SCDCalibratePanels2");
60 "Workspace of Indexed Peaks");
63 auto mustBeNonNegative = std::make_shared<BoundedValidator<double>>();
64 mustBeNonNegative->setLower(0.0);
65 declareProperty(
"RecalculateUB",
true,
"Recalculate UB matrix using given lattice constants");
66 declareProperty(
"a",
EMPTY_DBL(), mustBeNonNegative,
67 "Lattice Parameter a (Leave empty to use lattice constants "
68 "in peaks workspace)");
69 declareProperty(
"b",
EMPTY_DBL(), mustBeNonNegative,
70 "Lattice Parameter b (Leave empty to use lattice constants "
71 "in peaks workspace)");
72 declareProperty(
"c",
EMPTY_DBL(), mustBeNonNegative,
73 "Lattice Parameter c (Leave empty to use lattice constants "
74 "in peaks workspace)");
75 declareProperty(
"alpha",
EMPTY_DBL(), mustBeNonNegative,
76 "Lattice Parameter alpha in degrees (Leave empty to use "
77 "lattice constants in peaks workspace)");
78 declareProperty(
"beta",
EMPTY_DBL(), mustBeNonNegative,
79 "Lattice Parameter beta in degrees (Leave empty to use "
80 "lattice constants in peaks workspace)");
81 declareProperty(
"gamma",
EMPTY_DBL(), mustBeNonNegative,
82 "Lattice Parameter gamma in degrees (Leave empty to use "
83 "lattice constants in peaks workspace)");
84 const std::string LATTICE(
"Lattice Constants");
85 setPropertyGroup(
"RecalculateUB", LATTICE);
86 setPropertyGroup(
"a", LATTICE);
87 setPropertyGroup(
"b", LATTICE);
88 setPropertyGroup(
"c", LATTICE);
89 setPropertyGroup(
"alpha", LATTICE);
90 setPropertyGroup(
"beta", LATTICE);
91 setPropertyGroup(
"gamma", LATTICE);
92 setPropertySettings(
"a", std::make_unique<EnabledWhenProperty>(
"RecalculateUB",
IS_DEFAULT));
93 setPropertySettings(
"b", std::make_unique<EnabledWhenProperty>(
"RecalculateUB",
IS_DEFAULT));
94 setPropertySettings(
"c", std::make_unique<EnabledWhenProperty>(
"RecalculateUB",
IS_DEFAULT));
95 setPropertySettings(
"alpha", std::make_unique<EnabledWhenProperty>(
"RecalculateUB",
IS_DEFAULT));
96 setPropertySettings(
"beta", std::make_unique<EnabledWhenProperty>(
"RecalculateUB",
IS_DEFAULT));
97 setPropertySettings(
"gamma", std::make_unique<EnabledWhenProperty>(
"RecalculateUB",
IS_DEFAULT));
107 const std::string CALIBRATION(
"Calibration Options");
111 declareProperty(
"CalibrateL1",
true,
"Change the L1(source to sample) distance");
112 declareProperty(
"SearchRadiusL1", 0.1, mustBeNonNegative,
113 "Search radius of delta L1 in meters, which is used to constrain optimization search space"
114 "when calibrating L1");
116 setPropertySettings(
"SearchRadiusL1", std::make_unique<EnabledWhenProperty>(
"CalibrateL1",
IS_EQUAL_TO,
"1"));
118 setPropertyGroup(
"CalibrateL1", CALIBRATION);
119 setPropertyGroup(
"SearchRadiusL1", CALIBRATION);
123 declareProperty(
"CalibrateBanks",
false,
"Calibrate position and orientation of each bank.");
125 "SearchRadiusTransBank", 5e-2, mustBeNonNegative,
126 "This is the search radius (in meter) when calibrating component translations, used to constrain optimization"
127 "search space when calibration translation of banks");
128 declareProperty(
"SearchradiusRotXBank", 1.0, mustBeNonNegative,
129 "This is the search radius (in deg) when calibrating component reorientation, used to constrain "
130 "optimization search space");
131 declareProperty(
"SearchradiusRotYBank", 1.0, mustBeNonNegative,
132 "This is the search radius (in deg) when calibrating component reorientation, used to constrain "
133 "optimization search space");
134 declareProperty(
"SearchradiusRotZBank", 1.0, mustBeNonNegative,
135 "This is the search radius (in deg) when calibrating component reorientation, used to constrain "
136 "optimization search space");
137 declareProperty(
"CalibrateSize",
false,
"Calibrate detector size for each bank.");
138 declareProperty(
"SearchRadiusSize", 0.0, mustBeNonNegative,
139 "This is the search radius (unit less) of scale factor around at value 1.0 "
140 "when calibrating component size if it is a rectangualr detector.");
141 declareProperty(
"FixAspectRatio",
true,
142 "If true, the scaling factor for detector along X- and Y-axis "
143 "must be the same. Otherwise, the 2 scaling factors are free.");
144 declareProperty(
"BankName",
"",
145 "If given, only the specified bank/component will be calibrated."
146 "Otherwise, all banks will be calibrated.");
149 setPropertySettings(
"BankName", std::make_unique<EnabledWhenProperty>(
"CalibrateBanks",
IS_EQUAL_TO,
"1"));
150 setPropertySettings(
"SearchRadiusTransBank",
151 std::make_unique<EnabledWhenProperty>(
"CalibrateBanks",
IS_EQUAL_TO,
"1"));
152 setPropertySettings(
"SearchradiusRotXBank",
153 std::make_unique<EnabledWhenProperty>(
"CalibrateBanks",
IS_EQUAL_TO,
"1"));
154 setPropertySettings(
"SearchradiusRotYBank",
155 std::make_unique<EnabledWhenProperty>(
"CalibrateBanks",
IS_EQUAL_TO,
"1"));
156 setPropertySettings(
"SearchradiusRotZBank",
157 std::make_unique<EnabledWhenProperty>(
"CalibrateBanks",
IS_EQUAL_TO,
"1"));
158 setPropertySettings(
"CalibrateSize", std::make_unique<EnabledWhenProperty>(
"CalibrateBanks",
IS_EQUAL_TO,
"1"));
159 setPropertySettings(
"SearchRadiusSize", std::make_unique<EnabledWhenProperty>(
"CalibrateSize",
IS_EQUAL_TO,
"1"));
160 setPropertySettings(
"FixAspectRatio", std::make_unique<EnabledWhenProperty>(
"CalibrateSize",
IS_EQUAL_TO,
"1"));
162 setPropertyGroup(
"CalibrateBanks", CALIBRATION);
163 setPropertyGroup(
"BankName", CALIBRATION);
164 setPropertyGroup(
"SearchRadiusTransBank", CALIBRATION);
165 setPropertyGroup(
"SearchradiusRotXBank", CALIBRATION);
166 setPropertyGroup(
"SearchradiusRotYBank", CALIBRATION);
167 setPropertyGroup(
"SearchradiusRotZBank", CALIBRATION);
168 setPropertyGroup(
"CalibrateSize", CALIBRATION);
169 setPropertyGroup(
"SearchRadiusSize", CALIBRATION);
170 setPropertyGroup(
"FixAspectRatio", CALIBRATION);
175 declareProperty(
"CalibrateT0",
false,
"Calibrate the T0 (initial TOF)");
176 declareProperty(
"SearchRadiusT0", 10.0, mustBeNonNegative,
177 "Search radius of T0 (in ms), used to constrain optimization search space");
179 setPropertySettings(
"SearchRadiusT0", std::make_unique<EnabledWhenProperty>(
"CalibrateT0",
IS_EQUAL_TO,
"1"));
181 setPropertyGroup(
"CalibrateT0", CALIBRATION);
182 setPropertyGroup(
"SearchRadiusT0", CALIBRATION);
186 declareProperty(
"TuneSamplePosition",
false,
"Fine tunning sample position");
187 declareProperty(
"SearchRadiusSamplePos", 0.1, mustBeNonNegative,
188 "Search radius of sample position change (in meters), used to constrain optimization search space");
190 setPropertySettings(
"SearchRadiusSamplePos",
191 std::make_unique<EnabledWhenProperty>(
"TuneSamplePosition",
IS_EQUAL_TO,
"1"));
193 setPropertyGroup(
"TuneSamplePosition", CALIBRATION);
194 setPropertyGroup(
"SearchRadiusSamplePos", CALIBRATION);
198 "The workspace containing the calibration table.");
200 const std::vector<std::string> detcalExts{
".DetCal",
".Det_Cal"};
203 "Path to an ISAW-style .detcal file to save.");
206 "Path to an Mantid .xml description(for LoadParameterFile) file to "
210 "Path to an .csv file which contains the Calibration Table");
212 const std::string OUTPUT(
"Output");
213 setPropertyGroup(
"OutputWorkspace", OUTPUT);
214 setPropertyGroup(
"DetCalFilename", OUTPUT);
215 setPropertyGroup(
"XmlFilename", OUTPUT);
216 setPropertyGroup(
"CSVFilename", OUTPUT);
220 declareProperty(
"VerboseOutput",
false,
"Toggle of child algorithm console output.");
221 declareProperty(
"ProfileL1",
false,
"Perform profiling of objective function with given input for L1");
222 declareProperty(
"ProfileBanks",
false,
"Perform profiling of objective function with given input for Banks");
223 declareProperty(
"ProfileT0",
false,
"Perform profiling of objective function with given input for T0");
224 declareProperty(
"ProfileL1T0",
false,
"Perform profiling of objective function along L1 and T0");
226 const std::string ADVCNTRL(
"Advanced Option");
227 setPropertyGroup(
"VerboseOutput", ADVCNTRL);
228 setPropertyGroup(
"ProfileL1", ADVCNTRL);
229 setPropertyGroup(
"ProfileBanks", ADVCNTRL);
230 setPropertyGroup(
"ProfileT0", ADVCNTRL);
231 setPropertyGroup(
"ProfileL1T0", ADVCNTRL);
234 auto mustBePositive = std::make_shared<Kernel::BoundedValidator<int>>();
235 mustBePositive->setLower(0);
236 declareProperty(
"MaxFitIterations", 500, mustBePositive,
237 "Stop after this number of iterations if a good fit is not found");
246 std::map<std::string, std::string> issues;
259 (!pws->sample().hasOrientedLattice())) {
260 issues[
"RecalculateUB"] =
"Lattice constants are needed for peak "
261 "workspace without a UB mattrix";
266 throw std::runtime_error(
"calibrationTableColumnTypes and calibrationTableColumnTypes have different size.");
294 bool calibrateBanks =
getProperty(
"CalibrateBanks");
295 bool tuneSamplePos =
getProperty(
"TuneSamplePosition");
302 const std::string DetCalFilename =
getProperty(
"DetCalFilename");
303 const std::string XmlFilename =
getProperty(
"XmlFilename");
304 const std::string CSVFilename =
getProperty(
"CSVFilename");
308 double sizesearchradius =
getProperty(
"SearchRadiusSize");
309 bool fixdetxyratio =
getProperty(
"FixAspectRatio");
315 std::vector<std::pair<std::string, bool>> criteria{{
"BankName",
true}};
316 m_pws->sort(criteria);
348 g_log.
notice() <<
"** Calibrating L1 (moderator) as requested\n";
352 if (calibrateBanks) {
353 g_log.
notice() <<
"** Calibrating L2 and orientation (bank) as requested\n";
354 optimizeBanks(m_pws, pws_original, docalibsize, sizesearchradius, fixdetxyratio);
357 if (calibrateL1 && calibrateBanks) {
358 g_log.
notice() <<
"** Calibrating L1 (moderator) after bank adjusted\n";
376 if (calibrateT0 && !calibrateL1) {
381 g_log.
notice() <<
"** Calibrating T0 only as requested\n";
385 if (tuneSamplePos && !calibrateL1) {
386 g_log.
notice() <<
"** Tunning sample position only as requested\n";
390 if (calibrateT0 && tuneSamplePos && !calibrateL1) {
391 g_log.
warning() <<
"** You have chosen to calibrate T0 and sample position while ignoring"
392 <<
" L1, which means an iterative search outside this calibration is needed"
393 <<
" in order to find the minimum.\n";
397 g_log.
notice() <<
"-- Generate calibration table\n";
398 Instrument_sptr instCalibrated = std::const_pointer_cast<Geometry::Instrument>(m_pws->getInstrument());
403 if (!XmlFilename.empty()) {
407 if (!DetCalFilename.empty()) {
411 if (!CSVFilename.empty()) {
431 double original_L1 = std::abs(pws->getInstrument()->getSource()->getPos().Z());
434 bool tuneSamplepos =
getProperty(
"TuneSamplePosition");
440 auto objf = std::make_shared<SCDCalibratePanels2ObjFunc>();
442 std::vector<double> tofs =
captureTOF(pws_original);
443 objf->setPeakWorkspace(pws,
"moderator", tofs);
444 fitL1_alg->setProperty(
"Function", std::dynamic_pointer_cast<IFunction>(objf));
447 std::ostringstream tie_str;
448 tie_str <<
"DeltaX=0.0,DeltaY=0.0,"
449 <<
"RotX=0.0,RotY=0.0,RotZ=0.0";
450 if (!tuneSamplepos) {
451 tie_str <<
",DeltaSampleX=0.0,DeltaSampleY=0.0,DeltaSampleZ=0.0";
454 tie_str <<
",DeltaT0=" <<
m_T0;
456 std::ostringstream constraint_str;
458 r_L1 = std::abs(r_L1);
459 constraint_str << -r_L1 <<
"<DeltaZ<" << r_L1;
463 r_dT0 = std::abs(r_dT0);
464 constraint_str <<
"," << -r_dT0 <<
"<DeltaT0<" << r_dT0;
467 double r_dsp =
getProperty(
"SearchRadiusSamplePos");
468 r_dsp = std::abs(r_dsp);
469 constraint_str <<
"," << -r_dsp <<
"<DeltaSampleX<" << r_dsp
470 <<
"," << -r_dsp <<
"<DeltaSampleY<" << r_dsp
471 <<
"," << -r_dsp <<
"<DeltaSampleZ<" << r_dsp;
474 fitL1_alg->setProperty(
"Ties", tie_str.str());
475 fitL1_alg->setProperty(
"Constraints", constraint_str.str());
476 fitL1_alg->setProperty(
"InputWorkspace", l1ws);
477 fitL1_alg->setProperty(
"CreateOutput",
true);
478 fitL1_alg->setProperty(
"Output",
"fit");
479 fitL1_alg->executeAsChildAlg();
482 std::ostringstream calilog;
483 double chi2OverDOF = fitL1_alg->getProperty(
"OutputChi2overDoF");
486 double dL1_optimized = rst->getRef<
double>(
"Value", 2);
489 double dT0_optimized = rst->getRef<
double>(
"Value", 6);
495 double dsx_optimized = rst->getRef<
double>(
"Value", 7);
496 double dsy_optimized = rst->getRef<
double>(
"Value", 8);
497 double dsz_optimized = rst->getRef<
double>(
"Value", 9);
501 pws->getInstrument()->getSource()->getName(), pws);
502 m_T0 = dT0_optimized;
504 "sample-position", pws);
506 int npks = pws->getNumberPeaks();
507 calilog <<
"-- Fit L1 results using " << npks <<
" peaks:\n"
508 <<
" dL1: " << dL1_optimized <<
" \n"
509 <<
" L1 " << original_L1 <<
" -> " << -pws->getInstrument()->getSource()->getPos().Z() <<
" \n"
510 <<
" dT0 = " <<
m_T0 <<
" (ms)\n"
511 <<
" dSamplePos = (" << dsx_optimized <<
"," << dsy_optimized <<
"," << dsz_optimized <<
")\n"
512 <<
" chi2/DOF = " << chi2OverDOF <<
"\n";
526 const bool &docalibsize,
const double &sizesearchradius,
527 const bool &fixdetxyratio) {
530 for (
int i = 0; i < static_cast<int>(
m_BankNames.size()); ++i) {
533 const std::string bankname = *std::next(
m_BankNames.begin(), i);
534 const std::string pwsBankiName =
"_pws_" + bankname;
540 std::stringstream ss;
541 ss <<
"i = " << i <<
" m bank name = " << bankname;
543 ss <<
" ... True ...";
545 ss <<
" ... Stop ...";
556 std::vector<double> tofs =
captureTOF(pwsBanki_original);
560 int nBankPeaks = pwsBanki->getNumberPeaks();
563 std::ostringstream msg_npeakCheckFail;
564 msg_npeakCheckFail <<
"-- Bank " << bankname <<
" have only " << nBankPeaks <<
" (<" <<
MINIMUM_PEAKS_PER_BANK
565 <<
") Peaks, skipping\n";
576 auto objf = std::make_shared<SCDCalibratePanels2ObjFunc>();
577 objf->setPeakWorkspace(pwsBanki, bankname, tofs);
578 fitBank_alg->setProperty(
"Function", std::dynamic_pointer_cast<IFunction>(objf));
582 double searchRadiusRotX =
getProperty(
"SearchradiusRotXBank");
583 searchRadiusRotX = std::abs(searchRadiusRotX);
584 double searchRadiusRotY =
getProperty(
"SearchradiusRotYBank");
585 searchRadiusRotY = std::abs(searchRadiusRotY);
586 double searchRadiusRotZ =
getProperty(
"SearchradiusRotZBank");
587 searchRadiusRotZ = std::abs(searchRadiusRotZ);
589 double searchRadiusTran =
getProperty(
"SearchRadiusTransBank");
590 searchRadiusTran = std::abs(searchRadiusTran);
592 std::ostringstream tie_str;
593 tie_str <<
"DeltaSampleX=0.0,DeltaSampleY=0.0,DeltaSampleZ=0.0,"
594 <<
"DeltaT0=" <<
m_T0;
595 std::ostringstream constraint_str;
598 tie_str <<
",RotX=0.0";
600 constraint_str << -searchRadiusRotX <<
"<RotX<" << searchRadiusRotX <<
",";
604 tie_str <<
",RotY=0.0";
606 constraint_str << -searchRadiusRotY <<
"<RotY<" << searchRadiusRotY <<
",";
610 tie_str <<
",RotZ=0.0";
612 constraint_str << -searchRadiusRotZ <<
"<RotZ<" << searchRadiusRotZ <<
",";
616 tie_str <<
",DeltaX=0.0,DeltaY=0.0,DeltaZ=0.0";
618 constraint_str << -searchRadiusTran <<
"<DeltaX<" << searchRadiusTran <<
","
619 << -searchRadiusTran <<
"<DeltaY<" << searchRadiusTran <<
","
620 << -searchRadiusTran <<
"<DeltaZ<" << searchRadiusTran;
626 std::shared_ptr<const Geometry::RectangularDetector> rectDet =
627 std::dynamic_pointer_cast<const Geometry::RectangularDetector>(comp);
631 std::ostringstream scaleconstraints;
632 std::ostringstream scaleties;
633 if (rectDet && docalibsize) {
635 scaleconstraints << scales.first - sizesearchradius <<
" <=ScaleX<" << scales.first + sizesearchradius;
637 scaleties <<
"ScaleX=ScaleY";
639 scaleconstraints <<
"," << scales.second - sizesearchradius <<
" <=ScaleY<" << scales.second + sizesearchradius;
643 scaleties <<
"ScaleX=" << scales.first <<
", ScaleY=" << scales.second;
647 std::string fitconstraint{constraint_str.str()};
648 if (scaleconstraints.str() !=
"") {
649 if (fitconstraint ==
"")
650 fitconstraint += scaleconstraints.str();
652 fitconstraint +=
"," + scaleconstraints.str();
654 std::string fittie{tie_str.str()};
655 if (scaleties.str() !=
"") {
657 fittie += scaleties.str();
659 fittie +=
"," + scaleties.str();
662 g_log.
information(
"Fitting " + bankname +
": constraint = " + fitconstraint +
"\n\t tie = " + fittie);
666 fitBank_alg->setProperty(
"Ties", fittie);
667 if (fitconstraint !=
"")
668 fitBank_alg->setProperty(
"Constraints", fitconstraint);
669 fitBank_alg->setProperty(
"InputWorkspace", wsBankCali);
670 fitBank_alg->setProperty(
"CreateOutput",
true);
671 fitBank_alg->setProperty(
"Output",
"fit");
674 fitBank_alg->executeAsChildAlg();
677 double chi2OverDOF = fitBank_alg->getProperty(
"OutputChi2overDoF");
679 double dx = rstFitBank->getRef<
double>(
"Value", 0);
680 double dy = rstFitBank->getRef<
double>(
"Value", 1);
681 double dz = rstFitBank->getRef<
double>(
"Value", 2);
682 double drx = rstFitBank->getRef<
double>(
"Value", 3);
683 double dry = rstFitBank->getRef<
double>(
"Value", 4);
684 double drz = rstFitBank->getRef<
double>(
"Value", 5);
685 double scalex = rstFitBank->getRef<
double>(
"Value", 10);
686 double scaley = rstFitBank->getRef<
double>(
"Value", 11);
689 std::string bn = bankname;
690 std::ostringstream calilog;
691 if (pws->getInstrument()->getName().compare(
"CORELLI") == 0) {
692 bn.append(
"/sixteenpack");
695 if (rectDet && docalibsize) {
704 V3D dtrans(dx, dy, dz);
705 V3D drots(drx, dry, drz);
706 calilog <<
"-- Fit " << bn <<
" results using " << nBankPeaks <<
" peaks:\n"
707 <<
" d(x,y,z) = " << dtrans <<
"\n"
708 <<
" r(x,y,z) = " << drots <<
"\n"
709 <<
" scale(x, y) = " <<
scalex <<
", " <<
scaley <<
" chi2/DOF = " << chi2OverDOF <<
"\n";
744 auto objf = std::make_shared<SCDCalibratePanels2ObjFunc>();
746 std::vector<double> tofs =
captureTOF(pws_original);
747 objf->setPeakWorkspace(pws,
"none", tofs);
748 fitT0_alg->setProperty(
"Function", std::dynamic_pointer_cast<IFunction>(objf));
751 std::ostringstream tie_str;
752 tie_str <<
"DeltaX=0.0,DeltaY=0.0,DeltaZ=0.0,"
753 <<
"RotX=0.0,RotY=0.0,RotZ=0.0,"
754 <<
"DeltaSampleX=0.0,DeltaSampleY=0.0,DeltaSampleZ=0.0";
755 std::ostringstream constraint_str;
757 r_dT0 = std::abs(r_dT0);
758 constraint_str << -r_dT0 <<
"<DeltaT0<" << r_dT0;
761 fitT0_alg->setProperty(
"Ties", tie_str.str());
762 fitT0_alg->setProperty(
"Constraints", constraint_str.str());
763 fitT0_alg->setProperty(
"InputWorkspace", t0ws);
764 fitT0_alg->setProperty(
"CreateOutput",
true);
765 fitT0_alg->setProperty(
"Output",
"fit");
766 fitT0_alg->executeAsChildAlg();
769 std::ostringstream calilog;
770 double chi2OverDOF = fitT0_alg->getProperty(
"OutputChi2overDoF");
772 double dT0_optimized = rst->getRef<
double>(
"Value", 6);
775 m_T0 = dT0_optimized;
776 int npks = pws->getNumberPeaks();
778 calilog <<
"-- Fit T0 results using " << npks <<
" peaks:\n"
779 <<
" dT0 = " <<
m_T0 <<
" (ms)\n"
780 <<
" chi2/DOF = " << chi2OverDOF <<
"\n";
797 auto objf = std::make_shared<SCDCalibratePanels2ObjFunc>();
799 std::vector<double> tofs =
captureTOF(pws_original);
800 objf->setPeakWorkspace(pws,
"none", tofs);
801 fitSamplePos_alg->setProperty(
"Function", std::dynamic_pointer_cast<IFunction>(objf));
804 std::ostringstream tie_str;
805 tie_str <<
"DeltaX=0.0,DeltaY=0.0,DeltaZ=0.0,"
806 <<
"RotX=0.0,RotY=0.0,RotZ=0.0,"
807 <<
"DeltaT0=" <<
m_T0;
808 std::ostringstream constraint_str;
809 double r_dsp =
getProperty(
"SearchRadiusSamplePos");
810 r_dsp = std::abs(r_dsp);
811 constraint_str << -r_dsp <<
"<DeltaSampleX<" << r_dsp <<
"," << -r_dsp <<
"<DeltaSampleY<" << r_dsp <<
"," << -r_dsp
812 <<
"<DeltaSampleZ<" << r_dsp;
815 fitSamplePos_alg->setProperty(
"Ties", tie_str.str());
816 fitSamplePos_alg->setProperty(
"Constraints", constraint_str.str());
817 fitSamplePos_alg->setProperty(
"InputWorkspace", samplePosws);
818 fitSamplePos_alg->setProperty(
"CreateOutput",
true);
819 fitSamplePos_alg->setProperty(
"Output",
"fit");
820 fitSamplePos_alg->executeAsChildAlg();
823 std::ostringstream calilog;
824 double chi2OverDOF = fitSamplePos_alg->getProperty(
"OutputChi2overDoF");
826 double dsx_optimized = rst->getRef<
double>(
"Value", 7);
827 double dsy_optimized = rst->getRef<
double>(
"Value", 8);
828 double dsz_optimized = rst->getRef<
double>(
"Value", 9);
832 "sample-position", pws);
833 int npks = pws->getNumberPeaks();
835 calilog <<
"-- Tune SamplePos results using " << npks <<
" peaks:\n"
836 <<
" deltaSamplePos = (" << dsx_optimized <<
"," << dsy_optimized <<
"," << dsz_optimized <<
")\n"
837 <<
" chi2/DOF = " << chi2OverDOF <<
"\n";
861 (pws->sample().hasOrientedLattice())) {
881 calcUB_alg->setProperty(
"PeaksWorkspace", pws);
882 calcUB_alg->setProperty(
"a",
m_a);
883 calcUB_alg->setProperty(
"b",
m_b);
884 calcUB_alg->setProperty(
"c",
m_c);
885 calcUB_alg->setProperty(
"alpha",
m_alpha);
886 calcUB_alg->setProperty(
"beta",
m_beta);
887 calcUB_alg->setProperty(
"gamma",
m_gamma);
888 calcUB_alg->executeAsChildAlg();
893 idxpks_alg->setProperty(
"PeaksWorkspace", pws);
894 idxpks_alg->setProperty(
"RoundHKLs",
true);
895 idxpks_alg->setProperty(
"Tolerance", 0.15);
896 idxpks_alg->executeAsChildAlg();
908 fltpk_alg->setProperty(
"InputWorkspace", pws);
909 fltpk_alg->setProperty(
"FilterVariable",
"h^2+k^2+l^2");
910 fltpk_alg->setProperty(
"Operator",
">");
911 fltpk_alg->setProperty(
"FilterValue", 0.0);
912 fltpk_alg->setProperty(
"OutputWorkspace",
"pws_filtered");
913 fltpk_alg->executeAsChildAlg();
926 std::vector<double> tofs;
928 for (
int i = 0; i < pws->getNumberPeaks(); ++i) {
929 tofs.emplace_back(pws->getPeak(i).getTOF());
941 auto peaksWorkspace = std::dynamic_pointer_cast<DataObjects::PeaksWorkspace>(pws);
943 throw std::invalid_argument(
"a PeaksWorkspace is required to retrieve bank names");
944 int npeaks =
static_cast<int>(pws->getNumberPeaks());
945 for (
int i = 0; i < npeaks; ++i) {
946 std::string bname = peaksWorkspace->getPeak(i).getBankName();
961 const std::string &bankname,
962 const std::string &outputwsn) {
965 fltpk_alg->setProperty(
"InputWorkspace", pws);
966 fltpk_alg->setProperty(
"BankName", bankname);
967 fltpk_alg->setProperty(
"Criterion",
"=");
968 fltpk_alg->setProperty(
"OutputWorkspace", outputwsn);
969 fltpk_alg->executeAsChildAlg();
983 int npeaks = pws->getNumberPeaks();
991 auto &spectrum = mws->getSpectrum(0);
992 auto &xvector = spectrum.mutableX();
993 auto &yvector = spectrum.mutableY();
994 auto &evector = spectrum.mutableE();
997 double totalSigmaInt = 0.0;
998 for (
int i = 0; i < npeaks; ++i) {
999 totalSigmaInt += pws->getPeak(i).getSigmaIntensity();
1001 double totalInt = 0.0;
1002 for (
int i = 0; i < npeaks; ++i) {
1003 totalInt += pws->getPeak(i).getIntensity();
1005 double totalCnt = 0.0;
1006 for (
int i = 0; i < npeaks; ++i) {
1007 totalCnt += pws->getPeak(i).getBinCount();
1011 auto ubmatrix = pws->sample().getOrientedLattice().getUB();
1012 for (
int i = 0; i < npeaks; ++i) {
1014 V3D qv = ubmatrix * pws->getPeak(i).getIntHKL();
1018 if (totalSigmaInt > 0.0) {
1019 wgt = 1.0 / pws->getPeak(i).getSigmaIntensity();
1020 }
else if (totalInt > 0.0) {
1021 wgt = 1.0 / pws->getPeak(i).getIntensity();
1022 }
else if (totalCnt > 0.0) {
1023 wgt = 1.0 / pws->getPeak(i).getBinCount();
1026 for (
int j = 0; j < 3; ++j) {
1027 xvector[i * 3 + j] = i * 3 + j;
1028 yvector[i * 3 + j] = qv[j];
1029 evector[i * 3 + j] = wgt;
1051 double scalex,
double scaley,
const std::string &cmptName,
1057 mv_alg->setProperty(
"ComponentName", cmptName);
1058 mv_alg->setProperty(
"X", dx);
1059 mv_alg->setProperty(
"Y", dy);
1060 mv_alg->setProperty(
"Z", dz);
1061 mv_alg->setProperty(
"RelativePosition",
true);
1062 mv_alg->executeAsChildAlg();
1069 rot_alg->setProperty(
"ComponentName", cmptName);
1070 rot_alg->setProperty(
"X", 1.0);
1071 rot_alg->setProperty(
"Y", 0.0);
1072 rot_alg->setProperty(
"Z", 0.0);
1073 rot_alg->setProperty(
"Angle", drx);
1074 rot_alg->setProperty(
"RelativeRotation",
true);
1075 rot_alg->executeAsChildAlg();
1078 rot_alg->setProperty(
"ComponentName", cmptName);
1079 rot_alg->setProperty(
"X", 0.0);
1080 rot_alg->setProperty(
"Y", 1.0);
1081 rot_alg->setProperty(
"Z", 0.0);
1082 rot_alg->setProperty(
"Angle", dry);
1083 rot_alg->setProperty(
"RelativeRotation",
true);
1084 rot_alg->executeAsChildAlg();
1087 rot_alg->setProperty(
"ComponentName", cmptName);
1088 rot_alg->setProperty(
"X", 0.0);
1089 rot_alg->setProperty(
"Y", 0.0);
1090 rot_alg->setProperty(
"Z", 1.0);
1091 rot_alg->setProperty(
"Angle", drz);
1092 rot_alg->setProperty(
"RelativeRotation",
true);
1093 rot_alg->executeAsChildAlg();
1098 resizeAlg->initialize();
1099 resizeAlg->setProperty(
"Workspace", pws);
1100 resizeAlg->setProperty(
"ComponentName", cmptName);
1101 resizeAlg->setProperty(
"ScaleX",
scalex);
1102 resizeAlg->setProperty(
"ScaleY",
scaley);
1103 resizeAlg->execute();
1118 g_log.
notice() <<
"Generate a TableWorkspace to store calibration results.\n";
1127 V3D sourceRelPos = source->getRelativePos();
1131 sourceRow << instrument->getSource()->getName() << sourceRelPos.
X() << sourceRelPos.
Y() << sourceRelPos.
Z() << 1.0
1132 << 0.0 << 0.0 << 0.0 << 0.0 << 0.0;
1138 if (instrument->getName().compare(
"CORELLI") == 0)
1139 bankName.append(
"/sixteenpack");
1141 std::shared_ptr<const IComponent> bank = instrument->getComponentByName(bankName);
1143 Quat relRot = bank->getRelativeRot();
1144 V3D pos1 = bank->getRelativePos();
1147 double deg, xAxis, yAxis, zAxis;
1156 bankRow << bankName << pos1.
X() << pos1.
Y() << pos1.
Z() << xAxis << yAxis << zAxis << deg << scales.first
1183 const boost::container::flat_set<std::string> &AllBankNames,
1187 using boost::property_tree::ptree;
1192 parafile.put(
"<xmlattr>.instrument", instrument->getName());
1193 parafile.put(
"<xmlattr>.valid-from", instrument->getValidFromDate().toISO8601String());
1197 ptree src_dx, src_dy, src_dz;
1198 ptree src_dx_val, src_dy_val, src_dz_val;
1201 V3D sourceRelPos = source->getRelativePos();
1203 src_dx_val.put(
"<xmlattr>.val", sourceRelPos.
X());
1204 src_dy_val.put(
"<xmlattr>.val", sourceRelPos.
Y());
1205 src_dz_val.put(
"<xmlattr>.val", sourceRelPos.
Z());
1206 src_dx.put(
"<xmlattr>.name",
"x");
1207 src_dy.put(
"<xmlattr>.name",
"y");
1208 src_dz.put(
"<xmlattr>.name",
"z");
1209 src.put(
"<xmlattr>.name", source->getName());
1211 src_dx.add_child(
"value", src_dx_val);
1212 src_dy.add_child(
"value", src_dy_val);
1213 src_dz.add_child(
"value", src_dz_val);
1214 src.add_child(
"parameter", src_dx);
1215 src.add_child(
"parameter", src_dy);
1216 src.add_child(
"parameter", src_dz);
1218 parafile.add_child(
"component-link", src);
1223 ptree property_root;
1224 property_root.put(
"<xmlattr>.name", instrument->getName());
1227 tof0.put(
"<xmlattr>.name",
"T0");
1228 tof0_val.put(
"<xmlattr>.val",
m_T0);
1229 tof0.add_child(
"value", tof0_val);
1230 property_root.add_child(
"parameter", tof0);
1231 parafile.add_child(
"component-link", property_root);
1235 ptree samplePos_dx, samplePos_dy, samplePos_dz;
1236 ptree samplePos_dx_val, samplePos_dy_val, samplePos_dz_val;
1238 std::shared_ptr<const IComponent> sp = instrument->getComponentByName(
"sample-position");
1239 V3D sppos = sp->getRelativePos();
1240 samplePos_dx_val.put(
"<xmlattr>.val", sppos.
X());
1241 samplePos_dy_val.put(
"<xmlattr>.val", sppos.
Y());
1242 samplePos_dz_val.put(
"<xmlattr>.val", sppos.
Z());
1243 samplePos_dx.put(
"<xmlattr>.name",
"x");
1244 samplePos_dy.put(
"<xmlattr>.name",
"y");
1245 samplePos_dz.put(
"<xmlattr>.name",
"z");
1246 samplePos.put(
"<xmlattr>.name",
"sample-position");
1248 samplePos_dx.add_child(
"value", samplePos_dx_val);
1249 samplePos_dy.add_child(
"value", samplePos_dy_val);
1250 samplePos_dz.add_child(
"value", samplePos_dz_val);
1251 samplePos.add_child(
"parameter", samplePos_dx);
1252 samplePos.add_child(
"parameter", samplePos_dy);
1253 samplePos.add_child(
"parameter", samplePos_dz);
1255 parafile.add_child(
"component-link", samplePos);
1258 for (
auto bankName : AllBankNames) {
1260 if (instrument->getName().compare(
"CORELLI") == 0)
1261 bankName.append(
"/sixteenpack");
1263 std::shared_ptr<const IComponent> bank = instrument->getComponentByName(bankName);
1265 Quat relRot = bank->getRelativeRot();
1267 V3D pos1 = bank->getRelativePos();
1272 ptree bank_dx, bank_dy, bank_dz;
1273 ptree bank_dx_val, bank_dy_val, bank_dz_val;
1274 ptree bank_drotx, bank_droty, bank_drotz;
1275 ptree bank_drotx_val, bank_droty_val, bank_drotz_val;
1276 ptree bank_sx, bank_sy;
1277 ptree bank_sx_val, bank_sy_val;
1280 bank_dx_val.put(
"<xmlattr>.val", pos1.
X());
1281 bank_dy_val.put(
"<xmlattr>.val", pos1.
Y());
1282 bank_dz_val.put(
"<xmlattr>.val", pos1.
Z());
1283 bank_dx.put(
"<xmlattr>.name",
"x");
1284 bank_dy.put(
"<xmlattr>.name",
"y");
1285 bank_dz.put(
"<xmlattr>.name",
"z");
1287 bank_drotx_val.put(
"<xmlattr>.val", relRotAngles[0]);
1288 bank_droty_val.put(
"<xmlattr>.val", relRotAngles[1]);
1289 bank_drotz_val.put(
"<xmlattr>.val", relRotAngles[2]);
1290 bank_drotx.put(
"<xmlattr>.name",
"rotx");
1291 bank_droty.put(
"<xmlattr>.name",
"roty");
1292 bank_drotz.put(
"<xmlattr>.name",
"rotz");
1294 bank_sx_val.put(
"<xmlattr>.val", scales.first);
1295 bank_sy_val.put(
"<xmlattr>.val", scales.second);
1296 bank_sx.put(
"<xmlattr>.name",
"scalex");
1297 bank_sy.put(
"<xmlattr>.name",
"scaley");
1299 bank_root.put(
"<xmlattr>.name", bankName);
1302 bank_dx.add_child(
"value", bank_dx_val);
1303 bank_dy.add_child(
"value", bank_dy_val);
1304 bank_dz.add_child(
"value", bank_dz_val);
1306 bank_drotx.add_child(
"value", bank_drotx_val);
1307 bank_droty.add_child(
"value", bank_droty_val);
1308 bank_drotz.add_child(
"value", bank_drotz_val);
1310 bank_sx.add_child(
"value", bank_sx_val);
1311 bank_sy.add_child(
"value", bank_sy_val);
1313 bank_root.add_child(
"parameter", bank_drotx);
1315 bank_root.add_child(
"parameter", bank_droty);
1316 bank_root.add_child(
"parameter", bank_drotz);
1317 bank_root.add_child(
"parameter", bank_dx);
1318 bank_root.add_child(
"parameter", bank_dy);
1319 bank_root.add_child(
"parameter", bank_dz);
1320 bank_root.add_child(
"parameter", bank_sx);
1321 bank_root.add_child(
"parameter", bank_sy);
1323 parafile.add_child(
"component-link", bank_root);
1327 root.add_child(
"parameter-file", parafile);
1329 g_log.
notice() <<
"\tSaving parameter file as " << FileName <<
"\n";
1330 boost::property_tree::write_xml(FileName, root, std::locale(),
1331 boost::property_tree::xml_writer_settings<std::string>(
' ', 2));
1345 boost::container::flat_set<std::string> &AllBankName,
1346 std::shared_ptr<Instrument> &instrument,
double T0) {
1347 g_log.
notice() <<
"Saving DetCal file in " << filename <<
"\n";
1349 bool tuneSamplePos =
getProperty(
"TuneSamplePosition");
1350 if (tuneSamplePos) {
1352 <<
"DetCal format cannot retain sample position info, therefore the calibrated "
1353 <<
"sample position will be lost if DetCal format is the only output!\n";
1357 const size_t number_spectra = instrument->getNumberDetectors();
1360 wksp->setInstrument(instrument);
1361 wksp->rebuildSpectraMapping(
true );
1364 std::vector<std::string> banknames(AllBankName.begin(), AllBankName.end());
1368 alg->setProperty(
"InputWorkspace", wksp);
1369 alg->setProperty(
"Filename", filename);
1370 alg->setProperty(
"TimeOffset", T0);
1371 alg->setProperty(
"BankNames", banknames);
1372 alg->executeAsChildAlg();
1383 alg->setProperty(
"InputWorkspace", tws);
1384 alg->setProperty(
"Filename", FileName);
1385 alg->setPropertyValue(
"CommentIndicator",
"#");
1386 alg->setPropertyValue(
"Separator",
"CSV");
1387 alg->setProperty(
"ColumnHeader",
true);
1388 alg->setProperty(
"AppendToFile",
false);
1389 alg->executeAsChildAlg();
1400 g_log.
notice() <<
"START of profiling objective func along L1\n";
1410 std::ostringstream msgrst;
1411 msgrst.precision(12);
1412 msgrst <<
"dL1\tresidual\n";
1415 auto objf = std::make_shared<SCDCalibratePanels2ObjFunc>();
1417 std::vector<double> tofs =
captureTOF(pws_original);
1418 objf->setPeakWorkspace(pws,
"moderator", tofs);
1421 const int n_peaks = pws->getNumberPeaks();
1422 std::unique_ptr<double[]> target(
new double[n_peaks * 3]);
1425 auto ubmatrix = pws->sample().getOrientedLattice().getUB();
1426 for (
int i = 0; i < n_peaks; ++i) {
1427 V3D qv = ubmatrix * pws->getPeak(i).getIntHKL();
1429 for (
int j = 0; j < 3; ++j) {
1430 target[i * 3 + j] = qv[j];
1434 double xValues[7] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
1438 double deltaL1 = -4e-2;
1439 while (deltaL1 < 4e-2) {
1440 std::unique_ptr<double[]> out(
new double[n_peaks * 3]);
1441 objf->setParameter(
"DeltaZ", deltaL1);
1442 objf->setParameter(
"DeltaT0", 0.0);
1443 objf->function1D(out.get(), xValues, 1);
1446 double residual = 0.0;
1447 for (
int i = 0; i < n_peaks * 3; ++i) {
1448 residual += (out[i] - target[i]) * (out[i] - target[i]);
1450 residual = std::sqrt(residual) / (n_peaks - 1);
1452 msgrst << deltaL1 <<
"\t" << residual <<
"\n";
1455 g_log.
notice() << deltaL1 <<
" -- " << residual <<
"\n";
1463 auto filenamebase = boost::filesystem::temp_directory_path();
1464 filenamebase /= boost::filesystem::unique_path(
"profileSCDCalibratePanels2_L1.csv");
1465 std::ofstream profL1File;
1466 profL1File.open(filenamebase.string());
1467 profL1File << msgrst.str();
1470 << filenamebase <<
"\n"
1471 <<
"END of profiling objective func along L1\n";
1482 g_log.
notice() <<
"START of profiling all banks along six degree of freedom\n";
1493 for (
int bankIndex = 0; bankIndex < static_cast<int>(
m_BankNames.size()); ++bankIndex) {
1496 const std::string bankname = *std::next(
m_BankNames.begin(), bankIndex);
1497 const std::string pwsBankiName =
"_pws_" + bankname;
1503 std::vector<double> tofs =
captureTOF(pwsBanki_original);
1507 int nBankPeaks = pwsBanki->getNumberPeaks();
1510 std::ostringstream msg_npeakCheckFail;
1511 msg_npeakCheckFail <<
"-- Cannot profile Bank " << bankname <<
" have only " << nBankPeaks <<
" (<"
1519 std::ostringstream msgrst;
1520 msgrst.precision(12);
1521 msgrst <<
"dx\tdy\tdz\ttheta\tphi\trogang\tresidual\n";
1523 auto objf = std::make_shared<SCDCalibratePanels2ObjFunc>();
1524 objf->setPeakWorkspace(pwsBanki, bankname, tofs);
1526 const int n_peaks = pwsBanki->getNumberPeaks();
1527 std::unique_ptr<double[]> target(
new double[n_peaks * 3]);
1529 auto ubmatrix = pwsBanki->sample().getOrientedLattice().getUB();
1530 for (
int i = 0; i < n_peaks; ++i) {
1531 V3D qv = ubmatrix * pwsBanki->getPeak(i).getIntHKL();
1533 for (
int j = 0; j < 3; ++j) {
1534 target[i * 3 + j] = qv[j];
1538 double xValues[7] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
1541 for (
double dx = -1e-2; dx < 1e-2; dx += 2e-2 / 20.0) {
1543 for (
double dy = -1e-2; dy < 1e-2; dy += 2e-2 / 20.0) {
1545 for (
double dz = -1e-2; dz < 1e-2; dz += 2e-2 / 20.0) {
1547 for (
double theta = 0.0; theta <
PI; theta +=
PI / 20.0) {
1549 for (
double phi = 0.0; phi < 2 *
PI; phi += 2 *
PI / 20.0) {
1551 for (
double ang = -5.0; ang < 5.0; ang += 5.0 / 20.0) {
1554 std::unique_ptr<double[]> out(
new double[n_peaks * 3]);
1555 objf->setParameter(
"DeltaX", dx);
1556 objf->setParameter(
"DeltaY", dy);
1557 objf->setParameter(
"DeltaZ", dz);
1558 objf->setParameter(
"Theta", theta);
1559 objf->setParameter(
"Phi", phi);
1560 objf->setParameter(
"DeltaRotationAngle", ang);
1561 objf->setParameter(
"DeltaT0", 0.0);
1562 objf->function1D(out.get(), xValues, 1);
1564 double residual = 0.0;
1565 for (
int i = 0; i < n_peaks * 3; ++i) {
1566 residual += (out[i] - target[i]) * (out[i] - target[i]);
1568 residual = std::sqrt(residual) / (n_peaks - 6);
1570 msgrst << dx <<
"\t" << dy <<
"\t" << dz <<
"\t" << theta <<
"\t" << phi <<
"\t" << ang <<
"\t"
1571 << residual <<
"\n";
1574 g_log.
notice() <<
"--" << bankname <<
": " << residual <<
"\n";
1584 auto filenamebase = boost::filesystem::temp_directory_path();
1585 std::string fnbase =
"profileSCDCalibratePanels2_" + bankname +
".csv";
1586 filenamebase /= boost::filesystem::unique_path(fnbase);
1587 std::ofstream profBankFile;
1588 profBankFile.open(filenamebase.string());
1589 profBankFile << msgrst.str();
1590 profBankFile.close();
1593 std::ostringstream msg;
1594 msg <<
"Profile of " << bankname <<
" is saved at:\n"
1595 << filenamebase <<
"\n"
1596 <<
"END of profiling objective func for " << bankname <<
"\n";
1611 g_log.
notice() <<
"START of profiling objective func along T0\n";
1621 std::ostringstream msgrst;
1622 msgrst.precision(12);
1623 msgrst <<
"dT0\tresidual\n";
1626 auto objf = std::make_shared<SCDCalibratePanels2ObjFunc>();
1628 std::vector<double> tofs =
captureTOF(pws_original);
1629 objf->setPeakWorkspace(pws,
"none", tofs);
1632 const int n_peaks = pws->getNumberPeaks();
1633 std::unique_ptr<double[]> target(
new double[n_peaks * 3]);
1634 auto ubmatrix = pws->sample().getOrientedLattice().getUB();
1635 for (
int i = 0; i < n_peaks; ++i) {
1636 V3D qv = ubmatrix * pws->getPeak(i).getIntHKL();
1638 for (
int j = 0; j < 3; ++j) {
1639 target[i * 3 + j] = qv[j];
1643 double xValues[7] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
1646 double deltaT0 = -10;
1647 while (deltaT0 < 10) {
1648 std::unique_ptr<double[]> out(
new double[n_peaks * 3]);
1649 objf->setParameter(
"DeltaT0", deltaT0);
1650 objf->function1D(out.get(), xValues, 1);
1653 double residual = 0.0;
1654 for (
int i = 0; i < n_peaks * 3; ++i) {
1655 residual += (out[i] - target[i]) * (out[i] - target[i]);
1657 residual = std::sqrt(residual) / (n_peaks - 1);
1659 msgrst << deltaT0 <<
"\t" << residual <<
"\n";
1662 g_log.
notice() << deltaT0 <<
" -- " << residual <<
"\n";
1670 auto filenamebase = boost::filesystem::temp_directory_path();
1671 filenamebase /= boost::filesystem::unique_path(
"profileSCDCalibratePanels2_T0.csv");
1672 std::ofstream profL1File;
1673 profL1File.open(filenamebase.string());
1674 profL1File << msgrst.str();
1677 << filenamebase <<
"\n"
1678 <<
"END of profiling objective func along T0\n";
1689 g_log.
notice() <<
"START of profiling objective func along L1 and T0\n";
1695 g_log.
notice() <<
"deltaL1 -- deltaT0 -- residual\n";
1699 std::ostringstream msgrst;
1700 msgrst.precision(12);
1701 msgrst <<
"dL1\tdT0\tresidual\n";
1704 auto objf = std::make_shared<SCDCalibratePanels2ObjFunc>();
1706 std::vector<double> tofs =
captureTOF(pws_original);
1707 objf->setPeakWorkspace(pws,
"moderator", tofs);
1710 const int n_peaks = pws->getNumberPeaks();
1711 std::unique_ptr<double[]> target(
new double[n_peaks * 3]);
1712 auto ubmatrix = pws->sample().getOrientedLattice().getUB();
1713 for (
int i = 0; i < n_peaks; ++i) {
1714 V3D qv = ubmatrix * pws->getPeak(i).getIntHKL();
1716 for (
int j = 0; j < 3; ++j) {
1717 target[i * 3 + j] = qv[j];
1721 double xValues[7] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
1724 for (
double deltaL1 = -4e-2; deltaL1 < 4e-2; deltaL1 += 1e-4) {
1725 for (
double deltaT0 = -4.0; deltaT0 < 4.0; deltaT0 += 1e-2) {
1726 std::unique_ptr<double[]> out(
new double[n_peaks * 3]);
1727 objf->setParameter(
"DeltaZ", deltaL1);
1728 objf->setParameter(
"DeltaT0", deltaT0);
1729 objf->function1D(out.get(), xValues, 1);
1732 double residual = 0.0;
1733 for (
int i = 0; i < n_peaks * 3; ++i) {
1734 residual += (out[i] - target[i]) * (out[i] - target[i]);
1736 residual = std::sqrt(residual) / (n_peaks - 2);
1739 g_log.
notice() << deltaL1 <<
" -- " << deltaT0 <<
" -- " << residual <<
"\n";
1742 msgrst << deltaL1 <<
"\t" << deltaT0 <<
"\t" << residual <<
"\n";
1747 auto filenamebase = boost::filesystem::temp_directory_path();
1748 filenamebase /= boost::filesystem::unique_path(
"profileSCDCalibratePanels2_L1T0.csv");
1749 std::ofstream profL1File;
1750 profL1File.open(filenamebase.string());
1751 profL1File << msgrst.str();
1756 << filenamebase <<
"\n"
1757 <<
"END of profiling objective func along L1 and T0\n";
1772std::pair<double, double>
1774 const std::string &bankname,
1777 std::pair<double, double> scales{1.0, 1.0};
1781 std::shared_ptr<const Geometry::RectangularDetector> rectDet =
1782 std::dynamic_pointer_cast<const Geometry::RectangularDetector>(comp);
1786 auto scalexparams =
pmap.getDouble(rectDet->getName(),
"scalex");
1787 auto scaleyparams =
pmap.getDouble(rectDet->getName(),
"scaley");
1788 if (!scalexparams.empty())
1789 scales.first = scalexparams[0];
1790 if (!scaleyparams.empty())
1791 scales.second = scaleyparams[0];
#define DECLARE_ALGORITHM(classname)
#define PARALLEL_START_INTERRUPT_REGION
Begins a block to skip processing is the algorithm has been interupted Note the end of the block if n...
#define PARALLEL_END_INTERRUPT_REGION
Ends a block to skip processing is the algorithm has been interupted Note the start of the block if n...
#define PARALLEL_FOR_IF(condition)
Empty definitions - to enable set your complier to enable openMP.
#define PARALLEL_CHECK_INTERRUPT_REGION
Adds a check after a Parallel region to see if it was interupted.
std::string getPropertyValue(const std::string &name) const override
Get the value of a property as a string.
TypedValue getProperty(const std::string &name) const override
Get the value of a property.
virtual std::shared_ptr< Algorithm > createChildAlgorithm(const std::string &name, const double startProgress=-1., const double endProgress=-1., const bool enableLogging=true, const int &version=-1)
Create a Child Algorithm.
@ OptionalSave
to specify a file to write to but an empty string is
TableRow represents a row in a TableWorkspace.
A property class for workspaces.
SCDCalibratePanels2 : Using input peakworkspace with indexation results to calibrate each individual ...
Mantid::API::IPeaksWorkspace_sptr removeUnindexedPeaks(const Mantid::API::IPeaksWorkspace_sptr &pws)
Remove unindexed peaks from workspace.
void adjustComponent(double dx, double dy, double dz, double drx, double dry, double drz, double scalex, double scaley, const std::string &cmptName, Mantid::API::IPeaksWorkspace_sptr &pws)
Helper functions for adjusting components.
void saveIsawDetCal(const std::string &filename, boost::container::flat_set< std::string > &AllBankName, std::shared_ptr< Geometry::Instrument > &instrument, double T0)
Save to ISAW type det calibration output for backward compatiblity.
Mantid::API::ITableWorkspace_sptr generateCalibrationTable(std::shared_ptr< Geometry::Instrument > &instrument, const Geometry::ParameterMap &pmap)
Generate a Table workspace to store the calibration results.
void optimizeT0(Mantid::API::IPeaksWorkspace_sptr pws, Mantid::API::IPeaksWorkspace_sptr pws_original)
Private function for calibrating T0.
void optimizeL1(Mantid::API::IPeaksWorkspace_sptr pws, Mantid::API::IPeaksWorkspace_sptr pws_original)
Private function for calibrating L1.
double m_a
unique vars for a given instance of calibration
void optimizeSamplePos(Mantid::API::IPeaksWorkspace_sptr pws, Mantid::API::IPeaksWorkspace_sptr pws_original)
Private function for fine tunning sample position.
const std::vector< std::string > calibrationTableColumnTypes
void getBankNames(const Mantid::API::IPeaksWorkspace_sptr &pws)
Private function for getting names of banks to be calibrated.
void profileBanks(Mantid::API::IPeaksWorkspace_sptr &pws, const Mantid::API::IPeaksWorkspace_sptr &pws_original)
Profiling obj func along six degree of freedom, which can be very slow.
void saveCalibrationTable(const std::string &FileName, Mantid::API::ITableWorkspace_sptr &tws)
Save the calibration table to a CSV file.
boost::container::flat_set< std::string > m_BankNames
const int MINIMUM_PEAKS_PER_BANK
static constexpr double Tolerance
void updateUBMatrix(const Mantid::API::IPeaksWorkspace_sptr &pws)
Update the UB matrix.
void optimizeBanks(Mantid::API::IPeaksWorkspace_sptr pws, const Mantid::API::IPeaksWorkspace_sptr &pws_original, const bool &docalibsize, const double &sizesearchradius, const bool &fixdetxyratio)
Private function for calibrating banks.
void parseLatticeConstant(const Mantid::API::IPeaksWorkspace_sptr &pws)
Private function dedicated for parsing lattice constant.
const std::vector< std::string > calibrationTableColumnNames
void saveXmlFile(const std::string &FileName, const boost::container::flat_set< std::string > &AllBankNames, std::shared_ptr< Geometry::Instrument > &instrument, const Geometry::ParameterMap &pmap)
Save to xml file for Mantid to load by manual crafting.
void exec() override
Overwrites Algorithm method.
void profileT0(Mantid::API::IPeaksWorkspace_sptr &pws, Mantid::API::IPeaksWorkspace_sptr pws_original)
Profile obj func along T0 axis.
Mantid::API::IPeaksWorkspace_sptr selectPeaksByBankName(const Mantid::API::IPeaksWorkspace_sptr &pws, const std::string &bankname, const std::string &outputwsn)
Helper function for selecting peaks based on given bank name.
std::map< std::string, std::string > validateInputs() override
Private validator for inputs.
std::pair< double, double > getRectangularDetectorScaleFactors(std::shared_ptr< Geometry::Instrument > &instrument, const std::string &bankname, const Geometry::ParameterMap &pmap)
Retrieve "scalex" and "scaley" from a workspace's parameter map if the component is rectangular detec...
std::string mCalibBankName
Mantid::API::MatrixWorkspace_sptr getIdealQSampleAsHistogram1D(const Mantid::API::IPeaksWorkspace_sptr &pws)
Helper function that calculates the ideal qSample based on integer HKL.
void profileL1(Mantid::API::IPeaksWorkspace_sptr &pws, Mantid::API::IPeaksWorkspace_sptr pws_original)
Profile related functions.
void profileL1T0(Mantid::API::IPeaksWorkspace_sptr &pws, Mantid::API::IPeaksWorkspace_sptr pws_original)
Profile obj func along L1 and T0 axis.
std::vector< double > captureTOF(const Mantid::API::IPeaksWorkspace_sptr &pws)
Cache TOF equivalent to those measured from experiment.
Class to implement UB matrix.
double alpha() const
Get lattice parameter.
double a(int nd) const
Get lattice parameter a1-a3 as function of index (0-2)
double c() const
Get lattice parameter.
double beta() const
Get lattice parameter.
double b() const
Get lattice parameter.
double gamma() const
Get lattice parameter.
IPropertyManager * setProperty(const std::string &name, const T &value)
Templated method to set the value of a PropertyWithValue.
The Logger class is in charge of the publishing messages from the framework through various channels.
void notice(const std::string &msg)
Logs at notice level.
void warning(const std::string &msg)
Logs at warning level.
void information(const std::string &msg)
Logs at information level.
void getAngleAxis(double &_deg, double &_ax0, double &_ax1, double &ax2) const
Extracts the angle of roatation and axis.
std::vector< double > getEulerAngles(const std::string &convention) const
Calculate the Euler angles that are equivalent to this Quaternion.
static T & Instance()
Return a reference to the Singleton instance, creating it if it does not already exist Creation is do...
constexpr double X() const noexcept
Get x.
constexpr double Y() const noexcept
Get y.
constexpr double Z() const noexcept
Get z.
std::shared_ptr< IPeaksWorkspace > IPeaksWorkspace_sptr
shared pointer to Mantid::API::IPeaksWorkspace
std::shared_ptr< ITableWorkspace > ITableWorkspace_sptr
shared pointer to Mantid::API::ITableWorkspace
std::shared_ptr< Workspace > Workspace_sptr
shared pointer to Mantid::API::Workspace
std::shared_ptr< MatrixWorkspace > MatrixWorkspace_sptr
shared pointer to the matrix workspace base class
std::shared_ptr< Workspace2D > Workspace2D_sptr
shared pointer to Mantid::DataObjects::Workspace2D
std::unique_ptr< T > create(const P &parent, const IndexArg &indexArg, const HistArg &histArg)
This is the create() method that all the other create() methods call.
std::shared_ptr< const IComponent > IComponent_const_sptr
Typdef of a shared pointer to a const IComponent.
std::shared_ptr< Instrument > Instrument_sptr
Shared pointer to an instrument object.
std::enable_if< std::is_pointer< Arg >::value, bool >::type threadSafe(Arg workspace)
Thread-safety check Checks the workspace to ensure it is suitable for multithreaded access.
constexpr double EMPTY_DBL() noexcept
Returns what we consider an "empty" double within a property.
Generate a tableworkspace to store the calibration results.
adjust instrument component position and orientation
: detector size scale at y-direction
@ Input
An input workspace.
@ Output
An output workspace.