Skip to main content

astroimsim_data/
psf.rs

1use std::iter::Flatten;
2use std::path::PathBuf;
3use std::slice::Iter;
4use serde::{Deserialize, Serialize};
5use uvex_fitrs::{Fits, FitsData, FitsDataArray, Hdu, HeaderValue};
6use astroimsim_geometry::coordinate_system::{CoordinateSystem, Coordinates};
7use astroimsim_geometry::grid2d::GRID2D;
8use astroimsim_geometry::points::Point;
9use astroimsim_geometry::grid2d::InterpolationData;
10//use crate::point_source::{PointSource, SourceList};
11
12#[derive(Debug,Clone,Serialize)]
13pub struct PSF {
14    pub path: PathBuf,
15    pub data: Vec<Vec<f32>>,
16    pub x_pixels: usize,
17    pub y_pixels: usize,
18    pub center: Point,
19    pub size: (f64,f64),
20}
21
22pub struct DataFile {
23    pub description: String,
24    pub path: PathBuf,
25    pub x_pixels: usize,
26    pub y_pixels: usize,
27}
28
29
30pub enum Load{
31    FromKey(String),
32    FromValue(f64)
33}
34
35
36pub enum FITSType{
37    thirtytwo,
38    sixtyfour,
39}
40
41impl PSF {
42    pub fn load_file(file: PathBuf, center:(Load, Load), size:(Load, Load), x_num:usize, y_num:usize) -> PSF {
43        println!("Loading {:?} into a DataFrame ",file);
44        let fits = Fits::open(file.clone()).expect("Failed to open FITS file");
45        let primary_hdu= fits.iter().next().expect("Couldn't find primary HDU");
46        let (mut data,shape) = match primary_hdu.read_data() {
47            FitsData::FloatingPoint32(FitsDataArray { shape, data }) => (data,shape),
48            _ => panic!("Could not unpack PSF data")
49        }; //TODO add support for f64 etc
50
51        let normalization:f32 = data.iter().sum();
52        let data:Vec<f32> = data.iter().map(|x|x/normalization).collect();
53        println!("PSF has been normalized to {:?}",data.iter().sum::<f32>());
54
55        assert_eq!(shape[0], x_num,"Diva down! Tried to load a file with data of the wrong x size"); //check that the data is the expected size
56        assert_eq!(shape[1], y_num,"Diva down! Tried to load a file with data of the wrong y size");
57
58
59        let center_x:f64 = PSF::load(center.0, &primary_hdu);
60        let center_y:f64 = PSF::load(center.1, &primary_hdu);
61        let size_x:f64 = PSF::load(size.0, &primary_hdu);
62        let size_y:f64 = PSF::load(size.1, &primary_hdu);
63        let data = data.chunks(x_num).map(|i| i.to_vec()).collect();
64      //  println!("{:?}",data);
65        PSF {
66            path: file,
67            data,
68            x_pixels: x_num,
69            y_pixels: y_num,
70            center: Point::new(center_x,center_y,Coordinates::ABSOLUTE),
71            size: (size_x,size_y),
72        }
73    }
74    pub fn snap_to_grid(&self, grid: &GRID2D) -> usize{
75        let index = grid.snap(self.center.clone());
76        index
77    }
78
79
80    pub fn repack_data(flat_data: Vec<f32>) -> Vec<Vec<f32>>{
81        flat_data.chunks(64).map(|i| i.to_vec()).collect()
82    }
83
84    fn load(thingy:Load, header:&Hdu)-> f64{
85        match thingy{
86            Load::FromKey(Key) => {
87                match header.value(&Key).expect("failed to get key") {
88                    HeaderValue::RealFloatingNumber(value)=> *value,
89                    _ => panic!("could not unpack FITS header value")
90                }}
91            Load::FromValue(value) => value
92        }
93    }
94
95}
96
97
98use std::time::Instant;
99
100use astroimsim_spectra::spectral_response::SpectralResponseCurve;
101
102
103#[derive(Clone,Debug)]
104pub struct SpatialEffect {
105    pub label: &'static str,
106    pub grid: GRID2D,
107    pub data: Vec<Vec<f64>>, //data must be fractional - implement percent later?
108    pub fits_path: &'static str,
109}
110
111impl SpatialEffect{
112    pub fn new_empty(label:&'static str, grid:GRID2D,fits_path:&'static str)-> SpatialEffect{
113        SpatialEffect{label,grid,data:vec![],fits_path}
114    }
115
116    pub fn from_matrix(label:&'static str, grid:GRID2D,fits_path:&'static str,data:Vec<Vec<f64>>)-> SpatialEffect{
117        assert_eq!(grid.y_num,data.len());
118        for row in &data{ assert_eq!(grid.x_num, row.len()); }
119        SpatialEffect{label,grid,fits_path,data}
120    }
121
122    pub fn load_data(&mut self){
123        println!("Loading {:?} into {:?}",self.fits_path, self.label);
124        let fits = Fits::open(self.fits_path).expect("Failed to open FITS file");
125        let primary_hdu= fits.iter().next().expect("Couldn't find primary HDU");
126        let (mut data,shape) = match primary_hdu.read_data() {
127            FitsData::FloatingPoint64(FitsDataArray { shape, data }) => (data,shape),
128            FitsData::FloatingPoint32(FitsDataArray { shape, data }) => {
129                let data = data.iter().map(|x|*x as f64).collect();
130                (data,shape)},
131            _ => {panic!("huh? Couldn't load FITS file")}
132        };
133        assert_eq!(shape[0],self.grid.x_num,"FITS data had the wrong width");
134        assert_eq!(shape[1],self.grid.y_num,"FITS data had the wrong height");
135        let data:Vec<Vec<f64>> = data.chunks(self.grid.x_num).map(|v|v.to_vec()).collect();
136        self.data = data;
137
138    }
139    pub fn get_data_at_grid_index(&self, grid_number:usize)->f64{
140        let (x,y) = self.grid.xy_indices(grid_number);
141        self.data[y][x]
142    }
143
144    pub fn get_data(&self,point:&Point)->f64{
145        let interpolation_data = self.grid.interpolation_coefficients(point);
146        let sum:f64 = interpolation_data.corners
147            .iter()
148            .zip(interpolation_data.coefficients)
149            .map(|(point,coefficient)|{
150                self.get_data_at_grid_index(*point)*coefficient
151            }).sum();
152        sum/interpolation_data.normalization
153    }
154
155}
156
157
158
159
160
161
162
163
164
165