1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
|
// ************************************************************************************************
//
// BornAgain: simulate and fit reflection and scattering
//
//! @file Device/IO/ReadReflectometry.cpp
//! @brief Implements function readReflectometryTable.
//!
//! @homepage http://www.bornagainproject.org
//! @license GNU General Public License v3 or higher (see COPYING)
//! @copyright Forschungszentrum Jülich GmbH 2018
//! @authors Scientific Computing Group at MLZ (see CITATION, AUTHORS)
//
// ************************************************************************************************
#include "Device/IO/ReadReflectometry.h"
#include "Base/Axis/MakeScale.h"
#include "Base/Math/Numeric.h"
#include "Base/Util/Assert.h"
#include "Base/Util/StringUtil.h"
#include "Device/Data/Datafield.h"
#include "Device/IO/ImportSettings.h"
#include <algorithm>
#include <map>
#include <numbers>
using Base::String::split;
using Base::String::to_double;
using Base::String::to_int;
using Base::String::trim;
using std::numbers::pi;
Datafield Util::RW::readReflectometryTable(std::istream& s, const ImportSettings1D& p)
{
std::vector<int> ignorableLines = Base::String::expandNumberList(p.linesToSkip);
size_t maxCol = std::max({p.col_Q, p.col_R, p.col_dQ, p.col_sR, p.col_lambda});
std::string line;
int lineno = 0;
double2d_t rowsVec;
// Read numbers from file:
while (std::getline(s, line)) {
++lineno;
if (std::find(ignorableLines.begin(), ignorableLines.end(), lineno) != ignorableLines.end())
continue;
if (!p.headerPrefix.empty() && line.substr(0, p.headerPrefix.size()) == p.headerPrefix)
continue;
line = trim(line);
rowsVec.push_back(Base::String::parse_doubles(line)); // may throw
if (rowsVec.back().size() < maxCol)
throw std::runtime_error("Not enough entries in line " + std::to_string(lineno));
}
size_t arg_index = p.col_Q - 1; // p.col_Q >= 1
// sort
if (p.sort)
std::sort(rowsVec.begin(), rowsVec.end(),
[arg_index](const std::vector<double>& a, const std::vector<double>& b) {
return a[arg_index] < b[arg_index];
});
std::vector<double> QVec;
std::vector<double> RVec;
std::vector<double> sRVec;
double fac = 1.;
Coordinate outCoord = p.xCoord;
if (outCoord.unit() == "1/angstrom") {
fac = 10;
outCoord = {outCoord.name(), "1/nm"};
} else if (outCoord.unit() == "deg") {
fac = pi / 180.;
outCoord = {outCoord.name(), "rad"};
}
if (outCoord.name() == "2alpha") {
fac *= 0.5;
outCoord = {"alpha", outCoord.unit()};
}
for (const auto& row : rowsVec) {
const double arg = fac * row[arg_index];
if (p.rm_negative && arg < 0)
continue;
if (p.rm_duplications && !QVec.empty())
if (Numeric::almostEqual(QVec.back(), arg, 1))
continue;
QVec.push_back(arg);
RVec.push_back(row[p.col_R - 1]);
if (p.col_sR)
sRVec.push_back(row[p.col_sR - 1]);
}
return {std::vector<const Scale*>{newListScan(outCoord.label(), QVec)}, RVec, sRVec};
}
Datafield Util::RW::readMotofit(std::istream& s)
{
std::string line;
// First line
bool ok = (bool)std::getline(s, line);
if (!ok)
throw std::runtime_error("Empty input file");
if (trim(line) != "MFT")
throw std::runtime_error("File does not start with \"MFT\" line");
// Remaining header
int N = -1;
while (std::getline(s, line)) {
std::vector<std::string> splitted = split(line, ":");
if (splitted.size() < 2)
break; // end of header
if (splitted[0] == "Number of data points")
if (!to_int(trim(splitted[1]), &N))
throw std::runtime_error("Corrupt header line 'Number of data points'");
}
if (N == -1)
throw std::runtime_error("Missing header line 'Number of data points'");
// Blank line
if (!ok || !trim(line).empty())
throw std::runtime_error("Missing blank line after header");
// Tab header line
ok = (bool)std::getline(s, line);
if (!ok)
throw std::runtime_error("Missing tab header line");
if (line.size() < 80)
throw std::runtime_error("Tab header line shorter than expected");
if (trim(line.substr(0, 20)) != "q")
throw std::runtime_error("Unexpected entry 1 in tab header line");
if (trim(line.substr(20, 20)) != "refl")
throw std::runtime_error("Unexpected entry 2 in tab header line");
if (trim(line.substr(40, 20)) != "refl_err")
throw std::runtime_error("Unexpected entry 3 in tab header line");
if (trim(line.substr(60, 20)) != "q_res (FWHM)")
throw std::runtime_error("Unexpected entry 4 in tab header line");
// Data table
std::vector<double> Q(N), R(N), wR(N);
for (int i = 0; i < N; ++i) {
ok = (bool)std::getline(s, line);
if (!ok)
throw std::runtime_error("Missing line in data table");
if (line.size() < 80)
throw std::runtime_error("Line in data table shorter than expected");
ok = to_double(line.substr(0, 20), &(Q[i])) && to_double(line.substr(20, 20), &(R[i]))
&& to_double(line.substr(40, 20), &(wR[i]));
if (!ok)
throw std::runtime_error("Error in data line");
}
return {std::vector<const Scale*>{newListScan("q_z (1/nm)", Q)}, R, wR};
}
|
