1use std::collections::BTreeSet;
4use std::error::Error;
5use std::fmt::{Display, Formatter};
6
7use phasesmith_core::ConstantWavelengthInstrument;
8
9use crate::{
10 BackgroundError, DifferentiableBackground, LatticeBounds, ParameterBounds, ParameterError,
11 ParameterKey, ParameterSet, ParameterSpec, RietveldError, RietveldInput,
12 RietveldParameterError, RietveldStructuralLayout, RietveldStructuralSelection,
13 SamplePhysicsError,
14};
15
16#[derive(Clone, Copy, Debug, PartialEq, Eq, PartialOrd, Ord, Hash)]
18pub enum RietveldInstrumentParameter {
19 UDeg2,
21 VDeg2,
23 WDeg2,
25 XDeg,
27 YDeg,
29 WavelengthAngstrom,
31 ZeroShiftDeg,
33 SampleDisplacementMm,
35 DisplaceXMicrometre,
37 DisplaceYMicrometre,
39}
40
41impl RietveldInstrumentParameter {
42 #[must_use]
44 pub const fn as_str(self) -> &'static str {
45 match self {
46 Self::UDeg2 => "u_deg2",
47 Self::VDeg2 => "v_deg2",
48 Self::WDeg2 => "w_deg2",
49 Self::XDeg => "x_deg",
50 Self::YDeg => "y_deg",
51 Self::WavelengthAngstrom => "wavelength_angstrom",
52 Self::ZeroShiftDeg => "zero_shift_deg",
53 Self::SampleDisplacementMm => "sample_displacement_mm",
54 Self::DisplaceXMicrometre => "displace_x_micrometre",
55 Self::DisplaceYMicrometre => "displace_y_micrometre",
56 }
57 }
58}
59
60#[derive(Clone, Debug, Default, PartialEq, Eq)]
62pub struct RietveldParameterSelection {
63 pub structural: RietveldStructuralSelection,
65 pub instrument: Vec<RietveldInstrumentParameter>,
67 pub background: bool,
69 pub sample_physics: bool,
71}
72
73impl RietveldParameterSelection {
74 pub fn new(
81 structural: RietveldStructuralSelection,
82 instrument: Vec<RietveldInstrumentParameter>,
83 background: bool,
84 sample_physics: bool,
85 ) -> Result<Self, RietveldGeneralParameterError> {
86 let result = Self {
87 structural,
88 instrument,
89 background,
90 sample_physics,
91 };
92 result.validate()?;
93 Ok(result)
94 }
95
96 pub fn validate(&self) -> Result<(), RietveldGeneralParameterError> {
103 if self
104 .instrument
105 .iter()
106 .copied()
107 .collect::<BTreeSet<_>>()
108 .len()
109 != self.instrument.len()
110 {
111 return Err(RietveldGeneralParameterError::DuplicateInstrumentParameter);
112 }
113 Ok(())
114 }
115}
116
117#[derive(Clone, Debug, PartialEq)]
119pub struct RietveldParameterLayout {
120 parameters: ParameterSet,
121 structural: RietveldStructuralLayout,
122 structural_indices: Vec<usize>,
123 instrument: Vec<(RietveldInstrumentParameter, usize)>,
124 background_indices: Vec<usize>,
125 background_id: Option<String>,
126 sample_physics: Vec<SamplePhysicsMapping>,
127}
128
129#[derive(Clone, Debug, PartialEq, Eq)]
130struct SamplePhysicsMapping {
131 phase_index: usize,
132 name: String,
133 parameter_index: usize,
134}
135
136impl RietveldParameterLayout {
137 pub fn new(
144 input: &RietveldInput,
145 selection: &RietveldParameterSelection,
146 lattice_bounds: &[Option<LatticeBounds>],
147 ) -> Result<Self, RietveldGeneralParameterError> {
148 input.validate()?;
149 selection.validate()?;
150 if input.fixed_spectrum.is_some() {
151 if selection.structural.lattice {
152 return Err(RietveldGeneralParameterError::SpectrumLatticeRefinement);
153 }
154 if selection
155 .instrument
156 .contains(&RietveldInstrumentParameter::WavelengthAngstrom)
157 {
158 return Err(RietveldGeneralParameterError::SpectrumWavelengthRefinement);
159 }
160 }
161 let structural =
162 RietveldStructuralLayout::new(&input.phases, selection.structural, lattice_bounds)?;
163 let mut specs = Vec::new();
164 let mut instrument = Vec::new();
165 for selected in &selection.instrument {
166 let value = instrument_value(input, *selected)?;
167 let (unit, bounds, floor) = instrument_metadata(*selected)?;
168 let index = specs.len();
169 specs.push(ParameterSpec::new(
170 ParameterKey::new("instrument", "cw", selected.as_str())?,
171 value,
172 unit,
173 bounds,
174 value.abs().max(floor),
175 true,
176 )?);
177 instrument.push((*selected, index));
178 }
179 let mut background_indices = Vec::new();
180 let mut background_id = None;
181 if selection.background {
182 let background = input
183 .background
184 .as_ref()
185 .ok_or(RietveldGeneralParameterError::MissingBackground)?;
186 let names = background.parameter_names();
187 let coefficients = background.coefficients();
188 let bounds = background.parameter_bounds();
189 background_id = Some(background.background_id().to_owned());
190 for ((name, value), bounds) in names.into_iter().zip(coefficients).zip(bounds) {
191 let index = specs.len();
192 specs.push(ParameterSpec::new(
193 ParameterKey::new("background", background.background_id(), name)?,
194 value,
195 "intensity",
196 bounds,
197 value.abs().max(1.0),
198 true,
199 )?);
200 background_indices.push(index);
201 }
202 }
203 let mut sample_physics = Vec::new();
204 if selection.sample_physics {
205 for (phase_index, phase) in input.phases.iter().enumerate() {
206 let Some(model) = phase.sample_physics() else {
207 continue;
208 };
209 for parameter in model.parameters()? {
210 let parameter_index = specs.len();
211 specs.push(ParameterSpec::new(
212 ParameterKey::new("sample", phase.phase_id().as_str(), ¶meter.name)?,
213 parameter.value,
214 parameter.unit,
215 parameter.bounds,
216 parameter.scale,
217 true,
218 )?);
219 sample_physics.push(SamplePhysicsMapping {
220 phase_index,
221 name: parameter.name,
222 parameter_index,
223 });
224 }
225 }
226 }
227 let structural_indices = structural
228 .parameters()
229 .specs()
230 .iter()
231 .map(|spec| {
232 let index = specs.len();
233 specs.push(spec.clone());
234 index
235 })
236 .collect();
237 Ok(Self {
238 parameters: ParameterSet::new(specs)?,
239 structural,
240 structural_indices,
241 instrument,
242 background_indices,
243 background_id,
244 sample_physics,
245 })
246 }
247
248 #[must_use]
250 pub const fn parameters(&self) -> &ParameterSet {
251 &self.parameters
252 }
253
254 #[must_use]
256 pub const fn structural_layout(&self) -> &RietveldStructuralLayout {
257 &self.structural
258 }
259
260 pub fn structural_direction(
267 &self,
268 direction: &[f64],
269 ) -> Result<Vec<f64>, RietveldGeneralParameterError> {
270 self.validate_length(direction)?;
271 Ok(self
272 .structural_indices
273 .iter()
274 .map(|index| direction[*index])
275 .collect())
276 }
277
278 pub fn expand_structural_gradient(
285 &self,
286 structural: &[f64],
287 ) -> Result<Vec<f64>, RietveldGeneralParameterError> {
288 if structural.len() != self.structural_indices.len() {
289 return Err(RietveldGeneralParameterError::StructuralGradientLengthMismatch);
290 }
291 let mut result = vec![0.0; self.parameters.specs().len()];
292 for (value, index) in structural.iter().zip(&self.structural_indices) {
293 result[*index] = *value;
294 }
295 Ok(result)
296 }
297
298 pub(crate) fn instrument_indices(&self) -> &[(RietveldInstrumentParameter, usize)] {
299 &self.instrument
300 }
301
302 pub(crate) fn background_indices(&self) -> &[usize] {
303 &self.background_indices
304 }
305
306 pub(crate) fn sample_physics_indices(&self) -> impl Iterator<Item = (usize, &str, usize)> {
307 self.sample_physics.iter().map(|mapping| {
308 (
309 mapping.phase_index,
310 mapping.name.as_str(),
311 mapping.parameter_index,
312 )
313 })
314 }
315
316 pub fn apply_values(
323 &self,
324 input: &RietveldInput,
325 values: &[f64],
326 ) -> Result<RietveldInput, RietveldGeneralParameterError> {
327 let current = self
328 .parameters
329 .specs()
330 .iter()
331 .map(ParameterSpec::value)
332 .collect::<Vec<_>>();
333 self.apply_value_change(input, ¤t, values)
334 }
335
336 pub fn apply_value_change(
347 &self,
348 input: &RietveldInput,
349 current_values: &[f64],
350 values: &[f64],
351 ) -> Result<RietveldInput, RietveldGeneralParameterError> {
352 self.validate_length(current_values)?;
353 self.validate_length(values)?;
354 for (spec, value) in self.parameters.specs().iter().zip(values) {
355 if !value.is_finite() || !spec.bounds().contains(*value) {
356 return Err(RietveldGeneralParameterError::Parameter(
357 ParameterError::ValueOutsideBounds {
358 key: spec.key().clone(),
359 value: *value,
360 },
361 ));
362 }
363 }
364 let structural_values = self
365 .structural_indices
366 .iter()
367 .map(|index| values[*index])
368 .collect::<Vec<_>>();
369 let current_structural_values = self
370 .structural_indices
371 .iter()
372 .map(|index| current_values[*index])
373 .collect::<Vec<_>>();
374 let mut updated = input.clone();
375 updated.phases = self.structural.apply_value_change(
376 &input.phases,
377 ¤t_structural_values,
378 &structural_values,
379 )?;
380 let mut wavelength = None;
381 for (parameter, index) in &self.instrument {
382 install_instrument_value(&mut updated, *parameter, values[*index])?;
383 if *parameter == RietveldInstrumentParameter::WavelengthAngstrom {
384 wavelength = Some(values[*index]);
385 }
386 }
387 if let Some(wavelength) = wavelength {
388 updated.phases = updated
389 .phases
390 .iter()
391 .map(|phase| phase.with_wavelength(wavelength))
392 .collect::<Result<Vec<_>, _>>()?;
393 }
394 if !self.background_indices.is_empty() {
395 let background = updated
396 .background
397 .as_ref()
398 .ok_or(RietveldGeneralParameterError::MissingBackground)?;
399 if Some(background.background_id()) != self.background_id.as_deref() {
400 return Err(RietveldGeneralParameterError::BackgroundIdentityMismatch);
401 }
402 let coefficients = self
403 .background_indices
404 .iter()
405 .map(|index| values[*index])
406 .collect::<Vec<_>>();
407 updated.background = Some(background.replace_coefficients(&coefficients)?);
408 }
409 for phase_index in 0..updated.phases.len() {
410 let mappings = self
411 .sample_physics
412 .iter()
413 .filter(|mapping| mapping.phase_index == phase_index)
414 .collect::<Vec<_>>();
415 if mappings.is_empty() {
416 continue;
417 }
418 let model = updated.phases[phase_index]
419 .sample_physics()
420 .ok_or_else(|| RietveldGeneralParameterError::MissingSamplePhysics {
421 phase_id: updated.phases[phase_index].phase_id().to_string(),
422 })?;
423 let replacements = mappings
424 .into_iter()
425 .map(|mapping| (mapping.name.clone(), values[mapping.parameter_index]))
426 .collect();
427 updated.phases[phase_index] = updated.phases[phase_index]
428 .replace_sample_physics(model.replace_parameters(&replacements)?);
429 }
430 updated.validate()?;
431 Ok(updated)
432 }
433
434 fn validate_length(&self, values: &[f64]) -> Result<(), RietveldGeneralParameterError> {
435 if values.len() != self.parameters.specs().len() {
436 return Err(RietveldGeneralParameterError::ValueLengthMismatch);
437 }
438 Ok(())
439 }
440}
441
442fn instrument_value(
443 input: &RietveldInput,
444 parameter: RietveldInstrumentParameter,
445) -> Result<f64, RietveldGeneralParameterError> {
446 Ok(match parameter {
447 RietveldInstrumentParameter::UDeg2 => input.instrument.u_deg2,
448 RietveldInstrumentParameter::VDeg2 => input.instrument.v_deg2,
449 RietveldInstrumentParameter::WDeg2 => input.instrument.w_deg2,
450 RietveldInstrumentParameter::XDeg => input.instrument.x_deg,
451 RietveldInstrumentParameter::YDeg => input.instrument.y_deg,
452 RietveldInstrumentParameter::WavelengthAngstrom => input.instrument.wavelength_angstrom,
453 RietveldInstrumentParameter::ZeroShiftDeg => input.position_correction.zero_shift_deg,
454 RietveldInstrumentParameter::SampleDisplacementMm => input
455 .position_correction
456 .bragg_brentano_mm
457 .map(|value| value.0)
458 .ok_or(RietveldGeneralParameterError::InstrumentGeometryMismatch)?,
459 RietveldInstrumentParameter::DisplaceXMicrometre => input
460 .position_correction
461 .debye_scherrer_micrometre
462 .map(|value| value.0)
463 .ok_or(RietveldGeneralParameterError::InstrumentGeometryMismatch)?,
464 RietveldInstrumentParameter::DisplaceYMicrometre => input
465 .position_correction
466 .debye_scherrer_micrometre
467 .map(|value| value.1)
468 .ok_or(RietveldGeneralParameterError::InstrumentGeometryMismatch)?,
469 })
470}
471
472fn instrument_metadata(
473 parameter: RietveldInstrumentParameter,
474) -> Result<(&'static str, ParameterBounds, f64), ParameterError> {
475 Ok(match parameter {
476 RietveldInstrumentParameter::WavelengthAngstrom => (
477 "angstrom",
478 ParameterBounds::new(f64::MIN_POSITIVE, f64::INFINITY)?,
479 0.1,
480 ),
481 RietveldInstrumentParameter::SampleDisplacementMm => {
482 ("millimetre", ParameterBounds::default(), 1.0e-2)
483 }
484 RietveldInstrumentParameter::DisplaceXMicrometre
485 | RietveldInstrumentParameter::DisplaceYMicrometre => {
486 ("micrometre", ParameterBounds::default(), 1.0e3)
487 }
488 RietveldInstrumentParameter::UDeg2
489 | RietveldInstrumentParameter::VDeg2
490 | RietveldInstrumentParameter::WDeg2 => ("degree^2", ParameterBounds::default(), 1.0e-4),
491 RietveldInstrumentParameter::XDeg | RietveldInstrumentParameter::YDeg => {
492 ("degree", ParameterBounds::default(), 1.0e-3)
493 }
494 RietveldInstrumentParameter::ZeroShiftDeg => ("degree", ParameterBounds::default(), 5.0e-2),
495 })
496}
497
498fn install_instrument_value(
499 input: &mut RietveldInput,
500 parameter: RietveldInstrumentParameter,
501 value: f64,
502) -> Result<(), RietveldGeneralParameterError> {
503 let ConstantWavelengthInstrument {
504 wavelength_angstrom,
505 u_deg2,
506 v_deg2,
507 w_deg2,
508 x_deg,
509 y_deg,
510 } = &mut input.instrument;
511 match parameter {
512 RietveldInstrumentParameter::UDeg2 => *u_deg2 = value,
513 RietveldInstrumentParameter::VDeg2 => *v_deg2 = value,
514 RietveldInstrumentParameter::WDeg2 => *w_deg2 = value,
515 RietveldInstrumentParameter::XDeg => *x_deg = value,
516 RietveldInstrumentParameter::YDeg => *y_deg = value,
517 RietveldInstrumentParameter::WavelengthAngstrom => *wavelength_angstrom = value,
518 RietveldInstrumentParameter::ZeroShiftDeg => {
519 input.position_correction.zero_shift_deg = value;
520 }
521 RietveldInstrumentParameter::SampleDisplacementMm => {
522 let (_, radius) = input
523 .position_correction
524 .bragg_brentano_mm
525 .ok_or(RietveldGeneralParameterError::InstrumentGeometryMismatch)?;
526 input.position_correction.bragg_brentano_mm = Some((value, radius));
527 }
528 RietveldInstrumentParameter::DisplaceXMicrometre => {
529 let (_, y, radius) = input
530 .position_correction
531 .debye_scherrer_micrometre
532 .ok_or(RietveldGeneralParameterError::InstrumentGeometryMismatch)?;
533 input.position_correction.debye_scherrer_micrometre = Some((value, y, radius));
534 }
535 RietveldInstrumentParameter::DisplaceYMicrometre => {
536 let (x, _, radius) = input
537 .position_correction
538 .debye_scherrer_micrometre
539 .ok_or(RietveldGeneralParameterError::InstrumentGeometryMismatch)?;
540 input.position_correction.debye_scherrer_micrometre = Some((x, value, radius));
541 }
542 }
543 Ok(())
544}
545
546#[derive(Debug)]
548pub enum RietveldGeneralParameterError {
549 DuplicateInstrumentParameter,
551 MissingSamplePhysics {
553 phase_id: String,
555 },
556 MissingBackground,
558 InstrumentGeometryMismatch,
560 SpectrumWavelengthRefinement,
562 SpectrumLatticeRefinement,
564 BackgroundIdentityMismatch,
566 ValueLengthMismatch,
568 StructuralGradientLengthMismatch,
570 Parameter(ParameterError),
572 Structural(RietveldParameterError),
574 Rietveld(RietveldError),
576 Background(BackgroundError),
578 SamplePhysics(SamplePhysicsError),
580}
581
582impl Display for RietveldGeneralParameterError {
583 fn fmt(&self, formatter: &mut Formatter<'_>) -> std::fmt::Result {
584 match self {
585 Self::DuplicateInstrumentParameter => {
586 formatter.write_str("Rietveld instrument selections must be unique")
587 }
588 Self::MissingSamplePhysics { phase_id } => write!(
589 formatter,
590 "sample-physics model for phase {phase_id:?} changed under the parameter layout"
591 ),
592 Self::MissingBackground => {
593 formatter.write_str("background refinement requires an analytical background")
594 }
595 Self::InstrumentGeometryMismatch => formatter
596 .write_str("selected Rietveld position parameter does not match the geometry"),
597 Self::SpectrumWavelengthRefinement => {
598 formatter.write_str("fixed-spectrum component wavelengths cannot be refined")
599 }
600 Self::SpectrumLatticeRefinement => {
601 formatter.write_str("fixed-spectrum Rietveld lattice refinement is not supported")
602 }
603 Self::BackgroundIdentityMismatch => {
604 formatter.write_str("Rietveld background identity changed under the layout")
605 }
606 Self::ValueLengthMismatch => {
607 formatter.write_str("complete Rietveld value/direction length is wrong")
608 }
609 Self::StructuralGradientLengthMismatch => {
610 formatter.write_str("structural Rietveld gradient length is wrong")
611 }
612 Self::Parameter(error) => Display::fmt(error, formatter),
613 Self::Structural(error) => Display::fmt(error, formatter),
614 Self::Rietveld(error) => Display::fmt(error, formatter),
615 Self::Background(error) => Display::fmt(error, formatter),
616 Self::SamplePhysics(error) => Display::fmt(error, formatter),
617 }
618 }
619}
620
621impl Error for RietveldGeneralParameterError {
622 fn source(&self) -> Option<&(dyn Error + 'static)> {
623 match self {
624 Self::Parameter(error) => Some(error),
625 Self::Structural(error) => Some(error),
626 Self::Rietveld(error) => Some(error),
627 Self::Background(error) => Some(error),
628 Self::SamplePhysics(error) => Some(error),
629 Self::DuplicateInstrumentParameter
630 | Self::MissingSamplePhysics { .. }
631 | Self::MissingBackground
632 | Self::InstrumentGeometryMismatch
633 | Self::SpectrumWavelengthRefinement
634 | Self::SpectrumLatticeRefinement
635 | Self::BackgroundIdentityMismatch
636 | Self::ValueLengthMismatch
637 | Self::StructuralGradientLengthMismatch => None,
638 }
639 }
640}
641
642impl From<ParameterError> for RietveldGeneralParameterError {
643 fn from(value: ParameterError) -> Self {
644 Self::Parameter(value)
645 }
646}
647impl From<RietveldParameterError> for RietveldGeneralParameterError {
648 fn from(value: RietveldParameterError) -> Self {
649 Self::Structural(value)
650 }
651}
652impl From<RietveldError> for RietveldGeneralParameterError {
653 fn from(value: RietveldError) -> Self {
654 Self::Rietveld(value)
655 }
656}
657impl From<BackgroundError> for RietveldGeneralParameterError {
658 fn from(value: BackgroundError) -> Self {
659 Self::Background(value)
660 }
661}
662impl From<SamplePhysicsError> for RietveldGeneralParameterError {
663 fn from(value: SamplePhysicsError) -> Self {
664 Self::SamplePhysics(value)
665 }
666}