23#include <boost/algorithm/string/trim.hpp>
40 std::make_shared<InstrumentValidator>()),
41 "An input PeaksWorkspace with an instrument.");
43 declareProperty(
"AppendFile",
false,
44 "Append to file if true.\n"
45 "If false, new file (default).");
47 const std::vector<std::string> exts{
".peaks",
".integrate"};
48 declareProperty(std::make_unique<FileProperty>(
"Filename",
"",
FileProperty::Save, exts),
49 "Path to an ISAW-style peaks or integrate file to save.");
53 "An optional Workspace2D of profiles from integrating cylinder.");
55 declareProperty(
"RenumberPeaks",
false,
56 "If true, sequential peak numbers\n"
57 "If false, keep original numbering (default).");
64 std::string header =
"2 SEQN H K L COL ROW CHAN "
65 " L2 2_THETA AZ WL D IPK "
70 const auto &peaks = ws->getPeaks();
71 inst = ws->getInstrument();
73 throw std::runtime_error(
"No instrument in the Workspace. Cannot save DetCal file.");
74 const auto &detectorInfo = ws->detectorInfo();
78 using bankMap_t = std::map<int, std::vector<size_t>>;
79 using runMap_t = std::map<int, bankMap_t>;
80 std::set<int, std::less<int>> uniqueBanks;
83 std::string bankPart =
"bank";
84 if (
inst->getName() ==
"WISH")
85 bankPart =
"WISHpanel";
88 std::vector<IComponent_const_sptr> comps;
89 inst->getChildren(comps,
true);
91 for (
auto &comp : comps) {
92 std::string bankName = comp->getName();
93 boost::trim(bankName);
94 boost::erase_all(bankName, bankPart);
102 uniqueBanks.insert(bank);
105 for (
size_t i = 0; i < peaks.size(); ++i) {
106 const Peak &p = peaks[i];
112 if (bankName.size() <= 4) {
113 g_log.
information() <<
"Could not interpret bank number of peak " << i <<
"(" << bankName <<
")\n";
118 bankPart = bankName.substr(0, 4);
120 if (bankPart ==
"bank")
121 bankName = bankName.substr(4, bankName.size() - 4);
122 else if (bankPart ==
"WISHpanel")
123 bankName = bankName.substr(9, bankName.size() - 9);
127 runMap[run][bank].emplace_back(i);
130 header =
"2 SEQN H K L M N P COL ROW CHAN "
131 " L2 2_THETA AZ WL D IPK "
135 throw std::runtime_error(
"No instrument in PeaksWorkspace. Cannot save peaks file.");
137 if (bankPart !=
"bank" && bankPart !=
"WISHpanel" && bankPart !=
"?") {
138 std::ostringstream mess;
139 mess <<
"Detector module of type " << bankPart <<
" not supported in ISAWPeaks. Cannot save peaks file";
140 throw std::runtime_error(mess.str());
145 double beamline_norm;
147 inst->getInstrumentParameters(l1, beamline, beamline_norm, samplePos);
151 const bool renumber =
getProperty(
"RenumberPeaks");
154 if (!Poco::File(filename.c_str()).exists())
157 int appendPeakNumb = 0;
159 std::ifstream infile(filename.c_str());
161 while (!infile.eof())
163 getline(infile, line);
166 std::stringstream ss(line);
172 appendPeakNumb = std::max(peakNumber, appendPeakNumb);
177 out.open(filename.c_str(), std::ios::app);
178 appendPeakNumb = appendPeakNumb + 1;
180 out.open(filename.c_str());
182 const auto instrumentName =
inst->getName();
183 std::string facilityName;
188 <<
" not found at any defined facility. Setting facility "
190 facilityName =
"Unknown";
192 out <<
"Version: 2.0 Facility: " << facilityName;
193 out <<
" Instrument: " << instrumentName <<
" Date: ";
198 Types::Core::DateAndTime expDate =
inst->getValidFromDate() + 1.0;
199 out << expDate.toISO8601String();
204 out <<
"6 L1 T0_SHIFT\n";
205 out <<
"7 " << std::setw(10);
206 out << std::setprecision(4) << std::fixed << (l1 * 100);
207 out << std::setw(12) << std::setprecision(3) << std::fixed;
220 out <<
"4 DETNUM NROWS NCOLS WIDTH HEIGHT DEPTH DETD CenterX "
221 " CenterY CenterZ BaseX BaseY BaseZ UpX UpY "
224 for (
const auto bank : uniqueBanks) {
226 std::ostringstream mess;
227 if (bankPart ==
"bank")
228 mess <<
"bank" << bank;
229 else if (bankPart ==
"WISHpanel") {
230 mess <<
"WISHpanel" << std::setfill(
'0') << std::setw(2) << bank;
233 std::string bankName = mess.str();
235 std::shared_ptr<const IComponent> det =
inst->getComponentByName(bankName);
236 if (
inst->getName() ==
"CORELLI")
238 std::vector<Geometry::IComponent_const_sptr> children;
239 auto asmb = std::dynamic_pointer_cast<const Geometry::ICompAssembly>(
inst->getComponentByName(bankName));
240 asmb->getChildren(children,
false);
245 V3D center = det->getPos();
248 double detd = (center -
inst->getSample()->getPos()).norm();
251 sizeBanks(bankName, NCOLS, NROWS, xsize, ysize);
253 int midX = NCOLS / 2;
254 int midY = NROWS / 2;
263 out <<
"5 " << std::setw(6) << std::right << bank <<
" " << std::setw(6) << std::right << NROWS <<
" "
264 << std::setw(6) << std::right << NCOLS <<
" " << std::setw(7) << std::right << std::fixed
265 << std::setprecision(4) << 100.0 * xsize <<
" " << std::setw(7) << std::right << std::fixed
266 << std::setprecision(4) << 100.0 * ysize <<
" "
267 <<
" 0.2000 " << std::setw(6) << std::right << std::fixed << std::setprecision(2) << 100.0 * detd <<
" "
268 << std::setw(9) << std::right << std::fixed << std::setprecision(4) << 100.0 * center.
X() <<
" "
269 << std::setw(9) << std::right << std::fixed << std::setprecision(4) << 100.0 * center.
Y() <<
" "
270 << std::setw(9) << std::right << std::fixed << std::setprecision(4) << 100.0 * center.
Z() <<
" "
271 << std::setw(8) << std::right << std::fixed << std::setprecision(5) << base.
X() <<
" " << std::setw(8)
272 << std::right << std::fixed << std::setprecision(5) << base.
Y() <<
" " << std::setw(8) << std::right
273 << std::fixed << std::setprecision(5) << base.
Z() <<
" " << std::setw(8) << std::right << std::fixed
274 << std::setprecision(5) << up.
X() <<
" " << std::setw(8) << std::right << std::fixed << std::setprecision(5)
275 << up.
Y() <<
" " << std::setw(8) << std::right << std::fixed << std::setprecision(5) << up.
Z() <<
" \n";
278 g_log.
warning() <<
"Information about detector module " << bankName <<
" not found and recognised\n";
284 if (ws->getConvention() ==
"Crystallography")
289 int sequenceNumber = appendPeakNumb;
290 for (
const auto &runBankMap : runMap) {
292 const int run = runBankMap.first;
293 const auto &bankMap = runBankMap.second;
295 for (
const auto &bankIDs : bankMap) {
297 const int bank = bankIDs.first;
298 const auto &ids = bankIDs.second;
302 out <<
"0 NRUN DETNUM CHI PHI OMEGA MONCNT\n";
303 out <<
"1 " << std::setw(5) << run << std::setw(7) << std::right << bank;
307 Goniometer gon(peaks[ids[0]].getGoniometerMatrix());
310 double phi = angles[2];
311 double chi = angles[1];
312 double omega = angles[0];
314 out << std::setw(8) << std::fixed << std::setprecision(2) << chi <<
" ";
315 out << std::setw(8) << std::fixed << std::setprecision(2) << phi <<
" ";
316 out << std::setw(8) << std::fixed << std::setprecision(2) << omega <<
" ";
320 const size_t first_peak_index = ids[0];
321 const auto &first_peak = peaks[first_peak_index];
322 const double monct = first_peak.getMonitorCount();
323 out << std::setw(12) << static_cast<int>(monct) <<
'\n';
324 out << header <<
'\n';
327 for (
auto wi : ids) {
328 const auto &peak = peaks[wi];
331 std::string firstNumber =
"3";
336 out << firstNumber << std::setw(7) << sequenceNumber;
339 out << firstNumber << std::setw(7) << peak.getPeakNumber() + appendPeakNumb;
345 const V3D mod = peak.getIntMNP();
346 const auto intHKL = peak.getIntHKL();
353 out << std::setw(5) <<
Utils::round(qSign * peak.getH()) << std::setw(5)
358 out << std::setw(8) << std::fixed << std::setprecision(2) << static_cast<double>(peak.getCol()) <<
" ";
360 out << std::setw(8) << std::fixed << std::setprecision(2) << static_cast<double>(peak.getRow()) <<
" ";
362 out << std::setw(8) << std::fixed << std::setprecision(0) << peak.getTOF() <<
" ";
364 out << std::setw(9) << std::fixed << std::setprecision(3) << (peak.getL2() * 100.0) <<
" ";
367 const V3D dir = peak.getDetPos() -
inst->getSample()->getPos();
368 double scattering, azimuth;
372 scattering = dir.angle(
V3D(0.0, 0.0, 1.0));
377 azimuth = atan2(dir.Y(), dir.X());
379 out << std::setw(9) << std::fixed << std::setprecision(5) << scattering <<
" ";
381 out << std::setw(9) << std::fixed << std::setprecision(5) << azimuth <<
" ";
383 out << std::setw(10) << std::fixed << std::setprecision(6) << peak.getWavelength() <<
" ";
385 out << std::setw(9) << std::fixed << std::setprecision(4) << peak.getDSpacing() <<
" ";
387 out << std::setw(8) << std::fixed << std::setprecision(0) << int(peak.getBinCount()) <<
" ";
389 out << std::setw(10) << std::fixed << std::setprecision(2) << peak.getIntensity() <<
" ";
391 out << std::setw(7) << std::fixed << std::setprecision(2) << peak.getSigmaIntensity() <<
" ";
393 int thisReflag = 310;
394 out << std::setw(5) << thisReflag;
401 const auto &yValues = wsProfile2D->y(wi);
402 for (
size_t j = 0; j < yValues.size(); j++) {
403 out << std::setw(8) << static_cast<int>(yValues[j]);
404 if ((j + 1) % 10 == 0) {
406 if (j + 1 != yValues.size())
421 std::vector<Geometry::IComponent_const_sptr> children;
422 auto asmb = std::dynamic_pointer_cast<const Geometry::ICompAssembly>(parent);
423 asmb->getChildren(children,
false);
424 if (children[0]->
getName() ==
"sixteenpack") {
425 asmb = std::dynamic_pointer_cast<const Geometry::ICompAssembly>(children[0]);
427 asmb->getChildren(children,
false);
430 for (
const auto &col : children) {
431 auto asmb2 = std::dynamic_pointer_cast<const Geometry::ICompAssembly>(col);
432 std::vector<Geometry::IComponent_const_sptr> grandchildren;
433 asmb2->getChildren(grandchildren,
false);
434 for (
const auto &row : grandchildren) {
437 auto detID =
d->getID();
450 auto parent =
inst->getComponentByName(bankName);
451 if (parent->type() ==
"RectangularDetector") {
452 const auto RDet = std::dynamic_pointer_cast<const RectangularDetector>(parent);
453 const auto pixel = RDet->getAtXY(col, row);
454 return pixel->getPos();
456 std::vector<Geometry::IComponent_const_sptr> children;
457 auto asmb = std::dynamic_pointer_cast<const Geometry::ICompAssembly>(parent);
458 asmb->getChildren(children,
false);
459 if (children[0]->
getName() ==
"sixteenpack") {
460 asmb = std::dynamic_pointer_cast<const Geometry::ICompAssembly>(children[0]);
462 asmb->getChildren(children,
false);
466 if (
inst->getName() ==
"WISH")
467 col0 = (col % 2 == 0 ? col / 2 + 75 : (col - 1) / 2);
468 auto asmb2 = std::dynamic_pointer_cast<const Geometry::ICompAssembly>(children[col0]);
469 std::vector<Geometry::IComponent_const_sptr> grandchildren;
470 asmb2->getChildren(grandchildren,
false);
471 auto first = grandchildren[row - 1];
472 return first->getPos();
476 if (bankName ==
"None")
478 const auto parent =
inst->getComponentByName(bankName);
479 if (parent->type() ==
"RectangularDetector") {
480 const auto RDet = std::dynamic_pointer_cast<const RectangularDetector>(parent);
482 NCOLS = RDet->xpixels();
483 NROWS = RDet->ypixels();
484 xsize = RDet->xsize();
485 ysize = RDet->ysize();
487 std::vector<Geometry::IComponent_const_sptr> children;
488 auto asmb = std::dynamic_pointer_cast<const Geometry::ICompAssembly>(parent);
489 asmb->getChildren(children,
false);
490 if (children[0]->
getName() ==
"sixteenpack") {
491 asmb = std::dynamic_pointer_cast<const Geometry::ICompAssembly>(children[0]);
493 asmb->getChildren(children,
false);
495 const auto asmb2 = std::dynamic_pointer_cast<const Geometry::ICompAssembly>(children[0]);
496 std::vector<Geometry::IComponent_const_sptr> grandchildren;
497 asmb2->getChildren(grandchildren,
false);
498 NROWS =
static_cast<int>(grandchildren.size());
499 NCOLS =
static_cast<int>(children.size());
500 auto first = children[0];
501 auto last = children[NCOLS - 1];
502 xsize = first->getDistance(*last);
503 first = grandchildren[0];
504 last = grandchildren[NROWS - 1];
505 ysize = first->getDistance(*last);
509 for (
size_t i = 0; i < 3; i++) {
510 out << std::setw(12) << std::fixed << std::setprecision(6) << qSign * offset[i] <<
" ";
#define DECLARE_ALGORITHM(classname)
std::string getName(const IMDDimension &self)
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.
@ Save
to specify a file to write to, the file may or may not exist
bool hasProperty(const std::string &name) const
Does the property exist on the object.
HeldType getPropertyValueAsType(const std::string &name) const
Get the value of a property as the given TYPE.
This class stores information regarding an experimental run as a series of log entries.
A property class for workspaces.
Save a PeaksWorkspace to a ISAW-style ASCII .peaks file.
bool m_isModulatedStructure
Flag for writing modulated structures.
bool bankMasked(const Geometry::IComponent_const_sptr &parent, const Geometry::DetectorInfo &detectorInfo)
void writeOffsets(std::ofstream &out, double qSign, std::vector< double > offset)
void exec() override
Run the algorithm.
Kernel::V3D findPixelPos(const std::string &bankName, int col, int row)
find position for rectangular and non-rectangular
Geometry::Instrument_const_sptr inst
void sizeBanks(const std::string &bankName, int &NCOLS, int &NROWS, double &xsize, double &ysize)
Mantid::Kernel::V3D getIntMNP() const override
Return the int MNP vector.
int getRunNumber() const override
Return the run number this peak was measured at.
Structure describing a single-crystal peak.
std::string getBankName() const
Find the name of the bank that is the parent of the detector.
Geometry::DetectorInfo is an intermediate step towards a DetectorInfo that is part of Instrument-2....
bool isMasked(const size_t index) const
Returns true if the detector is masked.
size_t indexOf(const detid_t id) const
Returns the index of the detector with the given detector ID.
This class represents a detector - i.e.
Class to represent a particular goniometer setting, which is described by the rotation matrix.
std::vector< double > getEulerAngles(const std::string &convention="YZX")
Return Euler angles acording to a convention.
base class for Geometric IComponent
Exception for when an item is not found in a collection.
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.
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.
double normalize()
Make a normalized vector (return norm value)
constexpr double Y() const noexcept
Get y.
constexpr double Z() const noexcept
Get z.
std::shared_ptr< Workspace2D > Workspace2D_sptr
shared pointer to Mantid::DataObjects::Workspace2D
std::shared_ptr< PeaksWorkspace > PeaksWorkspace_sptr
Typedef for a shared pointer to a peaks workspace.
std::shared_ptr< const IComponent > IComponent_const_sptr
Typdef of a shared pointer to a const IComponent.
int convert(const std::string &A, T &out)
Convert a string into a number.
long round(double x)
Custom rounding method for a double->long because none is portable in C++ (!)
Peak indexing algorithm, which works by assigning multiple possible HKL values to each peak and then ...
@ Input
An input workspace.