BornAgain  1.19.0
Simulate and fit neutron and x-ray scattering at grazing incidence
RealSpace2DParacrystalUtils.cpp
Go to the documentation of this file.
1 // ************************************************************************************************
2 //
3 // BornAgain: simulate and fit reflection and scattering
4 //
5 //! @file GUI/coregui/Views/RealSpaceWidgets/RealSpace2DParacrystalUtils.cpp
6 //! @brief Defines RealSpaceBuilderUtils namespace
7 //!
8 //! @homepage http://www.bornagainproject.org
9 //! @license GNU General Public License v3 or higher (see COPYING)
10 //! @copyright Forschungszentrum Jülich GmbH 2018
11 //! @authors Scientific Computing Group at MLZ (see CITATION, AUTHORS)
12 //
13 // ************************************************************************************************
14 
18 
19 namespace {
20 void ResizeLatticePositions(std::vector<std::vector<double>>& lattice_positions, double l1,
21  double l2, double layer_size);
22 void FindLatticePositionsIndex(size_t& index, size_t& index_prev, int i, int j, int size,
23  double l_alpha);
24 std::pair<double, double> ComputePositionAlongPositiveLatticeVector(
25  const size_t index_prev, std::vector<std::vector<double>>& lattice_positions,
26  const IFTDistribution2D* pdf, double l, double l_xi, double l_alpha);
27 std::pair<double, double> ComputePositionAlongNegativeLatticeVector(
28  const size_t index_prev, std::vector<std::vector<double>>& lattice_positions,
29  const IFTDistribution2D* pdf, double l, double l_xi, double l_alpha);
30 std::pair<double, double>
31 ComputeLatticePosition(const size_t index_prev, int i, int j,
32  std::vector<std::vector<double>>& lattice_positions,
33  const IFTDistribution2D* pdf, double l, double l_xi, double l_alpha);
34 void ComputePositionsAlongLatticeVectorAxes(std::vector<std::vector<double>>& lattice_positions,
35  const IFTDistribution2D* pdf, double l, double l_xi,
36  double l_alpha);
37 void ComputePositionsInsideLatticeQuadrants(std::vector<std::vector<double>>& lattice_positions,
38  const IFTDistribution2D* pdf1,
39  const IFTDistribution2D* pdf2, double l1, double l2,
40  double l_xi, double l_alpha);
41 } // namespace
42 
44  const InterferenceFunction2DParaCrystal* p_iff, double layer_size)
45 {
46  auto& lattice = p_iff->lattice();
47  double l1 = lattice.length1();
48  double l2 = lattice.length2();
49  double alpha = lattice.latticeAngle();
50  double xi = lattice.rotationAngle();
51 
52  std::vector<std::vector<double>> lattice_positions;
53  ResizeLatticePositions(lattice_positions, l1, l2, layer_size);
54 
55  ComputePositionsAlongLatticeVectorAxes(lattice_positions, p_iff->pdf1(), l1, xi, 0);
56 
57  ComputePositionsAlongLatticeVectorAxes(lattice_positions, p_iff->pdf2(), l2, xi, alpha);
58 
59  ComputePositionsInsideLatticeQuadrants(lattice_positions, p_iff->pdf1(), p_iff->pdf2(), l1, l2,
60  xi, alpha);
61 
62  return lattice_positions;
63 }
64 
65 namespace {
66 void ResizeLatticePositions(std::vector<std::vector<double>>& lattice_positions, double l1,
67  double l2, double layer_size)
68 {
69  // Estimate the limit n1 and n2 of the integer multiple j and i of the lattice vectors l1 and l2
70  // required for populating particles correctly within the 3D model's boundaries
71  int n1 = 0, n2 = 0;
72  n1 = l1 == 0.0 ? 2 : static_cast<int>(layer_size * 2 / l1);
73  n2 = l2 == 0.0 ? 2 : static_cast<int>(layer_size * 2 / l2);
74 
75  n1 = std::max(n1, n2);
76 
77  lattice_positions.resize(static_cast<size_t>((2 * n1 + 1) * (2 * n1 + 1)));
78  for (auto& it : lattice_positions) {
79  it.resize(2);
80  }
81 
82  lattice_positions[0][0] = 0.0; // x coordinate of reference particle - at the origin
83  lattice_positions[0][1] = 0.0; // y coordinate of reference particle - at the origin
84 }
85 
86 void FindLatticePositionsIndex(size_t& index, size_t& index_prev, int i, int j, int size,
87  double l_alpha)
88 {
89  index = static_cast<size_t>(i * (2 * size + 1) + j);
90 
91  if (std::sin(l_alpha) == 0) // along l1
92  {
93  // particles along +l1 stored every odd iter (e.g. 1,3,5...) index of lattice_positions
94  // particles along -l1 stored every even iter (e.g. 2,4,6...) index of lattice_positions
95  index_prev = static_cast<size_t>(i * (2 * size + 1));
96  if (j - 2 > 0)
97  index_prev = index - 2;
98  } else // along l2
99  {
100  // particles along +l2/-l2 stored every (odd/even iter)*(2*n1+1) index of lattice_positions
101  index_prev = static_cast<size_t>(j);
102  if (i - 2 > 0)
103  index_prev = index - static_cast<size_t>(2 * (2 * size + 1));
104  }
105 }
106 
107 std::pair<double, double> ComputePositionAlongPositiveLatticeVector(
108  const size_t index_prev, std::vector<std::vector<double>>& lattice_positions,
109  const IFTDistribution2D* pdf, double l, double l_xi, double l_alpha)
110 {
111  double gamma_pdf = pdf->gamma();
112  std::pair<double, double> sampleXYpdf = pdf->createSampler()->randomSample();
113 
114  double offset_x_pdf = sampleXYpdf.first;
115  double offset_y_pdf = sampleXYpdf.second;
116 
117  double x = lattice_positions[index_prev][0] + l * std::cos(l_xi + l_alpha)
118  + offset_x_pdf * std::cos(gamma_pdf + l_xi)
119  + offset_y_pdf * std::cos(M_PI_2 + gamma_pdf + l_xi); // x coordinate
120  double y = lattice_positions[index_prev][1] + l * std::sin(l_xi + l_alpha)
121  + offset_x_pdf * std::sin(gamma_pdf + l_xi)
122  + offset_y_pdf * std::sin(M_PI_2 + gamma_pdf + l_xi); // y coordinate
123 
124  return std::make_pair(x, y);
125 }
126 
127 std::pair<double, double> ComputePositionAlongNegativeLatticeVector(
128  const size_t index_prev, std::vector<std::vector<double>>& lattice_positions,
129  const IFTDistribution2D* pdf, double l, double l_xi, double l_alpha)
130 {
131  double gamma_pdf = pdf->gamma();
132  std::pair<double, double> sampleXYpdf = pdf->createSampler()->randomSample();
133 
134  double offset_x_pdf = sampleXYpdf.first;
135  double offset_y_pdf = sampleXYpdf.second;
136 
137  double x = lattice_positions[index_prev][0] - l * std::cos(l_xi + l_alpha)
138  + offset_x_pdf * std::cos(gamma_pdf + l_xi)
139  + offset_y_pdf * std::cos(M_PI_2 + gamma_pdf + l_xi); // x coordinate
140  double y = lattice_positions[index_prev][1] - l * std::sin(l_xi + l_alpha)
141  + offset_x_pdf * std::sin(gamma_pdf + l_xi)
142  + offset_y_pdf * std::sin(M_PI_2 + gamma_pdf + l_xi); // y coordinate
143 
144  return std::make_pair(x, y);
145 }
146 
147 std::pair<double, double>
148 ComputeLatticePosition(const size_t index_prev, int i, int j,
149  std::vector<std::vector<double>>& lattice_positions,
150  const IFTDistribution2D* pdf, double l, double l_xi, double l_alpha)
151 {
152  if (std::sin(l_alpha) == 0) {
153  if (!(j % 2 == 0)) // along +l1
154  return ComputePositionAlongPositiveLatticeVector(index_prev, lattice_positions, pdf, l,
155  l_xi, 0);
156  else // along -l1
157  return ComputePositionAlongNegativeLatticeVector(index_prev, lattice_positions, pdf, l,
158  l_xi, 0);
159  } else {
160  if (!(i % 2 == 0)) // along +l2
161  return ComputePositionAlongPositiveLatticeVector(index_prev, lattice_positions, pdf, l,
162  l_xi, l_alpha);
163  else // along -l2
164  return ComputePositionAlongNegativeLatticeVector(index_prev, lattice_positions, pdf, l,
165  l_xi, l_alpha);
166  }
167 }
168 
169 void ComputePositionsAlongLatticeVectorAxes(std::vector<std::vector<double>>& lattice_positions,
170  const IFTDistribution2D* pdf, double l, double l_xi,
171  double l_alpha)
172 {
173  int n = static_cast<int>((std::sqrt(lattice_positions.size()) - 1) / 2);
174 
175  size_t index = 0; // lattice_positions index for current particle
176  size_t index_prev = 0; // lattice_positions index for previous particle
177 
178  std::pair<double, double> xy;
179 
180  for (int iter = 1; iter <= 2 * n; ++iter) {
181 
182  int iterl1, iterl2;
183 
184  if (std::sin(l_alpha) == 0) {
185  iterl1 = iter;
186  iterl2 = 0;
187  } else {
188  iterl1 = 0;
189  iterl2 = iter;
190  }
191 
192  // The 2*n+1 particles that are situated ONLY along the l1 axis (both +/- axes)
193  // are stored in i = 1,2,3,...2*n1 indices of lattice_positions
194 
195  // The 2*n+1 particles that are situated ONLY along the l2 axis (both +/- axes)
196  // are stored every i*(2*n1+1) index of lattice_positions
197 
198  FindLatticePositionsIndex(index, index_prev, iterl2, iterl1, n, l_alpha);
199  xy = ComputeLatticePosition(index_prev, iterl2, iterl1, lattice_positions, pdf, l, l_xi,
200  l_alpha);
201 
202  lattice_positions[index][0] = xy.first; // x coordinate
203  lattice_positions[index][1] = xy.second; // y coordinate
204  }
205 }
206 
207 void ComputePositionsInsideLatticeQuadrants(std::vector<std::vector<double>>& lattice_positions,
208  const IFTDistribution2D* pdf1,
209  const IFTDistribution2D* pdf2, double l1, double l2,
210  double l_xi, double l_alpha)
211 {
212  int n = static_cast<int>((std::sqrt(lattice_positions.size()) - 1) / 2);
213 
214  size_t index = 0; // lattice_positions index for current particle
215  size_t index_prev = 0; // lattice_positions index for previous particle
216 
217  std::pair<double, double> xy_l1, xy_l2;
218 
219  for (int i = 1; i <= 2 * n; ++i) {
220  for (int j = 1; j <= 2 * n; ++j) {
221  FindLatticePositionsIndex(index, index_prev, i, j, n, 0);
222  xy_l1 = ComputeLatticePosition(index_prev, i, j, lattice_positions, pdf1, l1, l_xi, 0);
223 
224  FindLatticePositionsIndex(index, index_prev, i, j, n, l_alpha);
225  xy_l2 = ComputeLatticePosition(index_prev, i, j, lattice_positions, pdf2, l2, l_xi,
226  l_alpha);
227 
228  lattice_positions[index][0] = (xy_l1.first + xy_l2.first) / 2;
229  lattice_positions[index][1] = (xy_l1.second + xy_l2.second) / 2;
230  }
231  }
232 }
233 } // namespace
#define M_PI_2
Definition: Constants.h:45
Defines class InterferenceFunction2DParaCrystal.
Defines RealSpaceBuilderUtils namespace.
Defines class RealSpaceCanvas.
Interface for two-dimensional distributions in Fourier space.
virtual std::unique_ptr< IDistribution2DSampler > createSampler() const =0
double gamma() const
Interference function of a 2D paracrystal.
virtual double length1() const =0
std::vector< std::vector< double > > Compute2DParacrystalLatticePositions(const InterferenceFunction2DParaCrystal *, double layer_size)