smlm_imp/
lib.rs

1extern crate smlm_sig_proc as sig_proc;
2extern crate smlm_tiff as tiff_wrap;
3
4#[cfg(feature = "plot")]
5extern crate gnuplot;
6
7pub use self::stats::Stats;
8pub use self::images::{OwnedImage, BorrowedImage};
9
10use tiff_wrap::{Tiff};
11
12pub use sig_proc::fft::c64;
13use sig_proc::fft::{fft_2d, fft_shift};
14
15#[cfg(feature = "plot")]
16use gnuplot::{Figure, AxesCommon, DataType};
17
18use std::ops::{AddAssign, Mul, Div};
19use std::path::Path;
20
21// mod convolution;
22pub mod circles;
23pub mod filter;
24pub mod images;
25pub mod stats;
26pub mod utils;
27
28pub trait Image
29{
30    type Data;
31    fn n_rows(&self) -> usize;
32    fn n_cols(&self) -> usize;
33    fn data(&self) -> &[Self::Data];
34
35    fn width(&self) -> usize
36    {
37        self.n_cols()
38    }
39
40    fn height(&self) -> usize
41    {
42        self.n_rows()
43    }
44
45    /// shape of image in rows x cols format
46    fn shape(&self) -> (usize, usize)
47    {
48        (self.n_rows(), self.n_cols())
49    }
50
51    /// n_rows * n_cols
52    fn size(&self) -> usize
53    {
54        self.n_rows() * self.n_cols()
55    }
56
57    /// get data with position
58    /// unsafe
59    fn get_at(&self, row: usize, col: usize) -> &Self::Data
60    {
61        let index = self.get_index(row, col);
62        self.get(index)
63    }
64
65    /// get data with index
66    /// unsafe
67    fn get(&self, idx: usize) -> &Self::Data
68    {
69        &self.data()[idx]
70    }
71
72    /// get data with indicies
73    /// unsafe
74    fn get_data(&self, indicies: &[usize]) -> Vec<&Self::Data>
75    {
76        indicies.iter().map(|x| self.get(*x)).collect()
77    }
78
79    /// get data index for row and column
80    fn get_index(&self, row: usize, col: usize) -> usize
81    {
82        utils::get_index(row, col, self.n_cols())
83    }
84
85    /// get coordinates for index
86    fn get_coords(&self, index: usize) -> (usize, usize)
87    {
88        utils::get_coords(index, self.n_cols())
89    }
90
91
92    //useful for things like frc
93    fn is_square(&self) -> bool
94    {
95        self.n_rows() == self.n_cols()
96    }
97
98    /// get data index for centre of the image
99    fn centre_index(&self) -> usize
100    {
101        let (row, col) = self.centre_position();
102        utils::get_index(row, col, self.n_cols())
103    }
104
105    /// get central row
106    fn centre_row(&self) -> usize
107    {
108        self.n_rows() / 2
109    }
110
111    /// get central column
112    fn centre_col(&self) -> usize
113    {
114        self.n_cols() / 2
115    }
116
117    /// get row and column for central pixel
118    fn centre_position(&self) -> (usize, usize)
119    {
120        (self.centre_row(), self.centre_col())
121    }
122
123    /// determine if image has a sqaure centre
124    /// ie 4 pixels make up central part
125    fn square_centre(&self) -> bool
126    {
127        self.n_rows() % 2 == 0 && self.n_rows() == self.n_cols()
128    }
129
130    ///indices for perimenter fo circle with specified radius using custom method
131    fn perimimeter_indicies_using<P: circles::Perimeter>(&self, radius: usize) -> Vec<usize>
132    where Self: Sized
133    {
134        P::get(self, radius as f64)
135    }
136
137    ///indices for perimenter fo circle with specified radius using known method
138    fn perimimeter_indicies_with(&self, method: circles::Method, radius: usize) -> Vec<usize>
139    where Self: Sized
140    {
141        method.get_perimeter(self, radius as f64)
142    }
143
144    ///indices for perimenter fo circle with specified radious
145    /// defaults to angle carve
146    fn perimimeter_indicies(&self, radius: usize) -> Vec<usize>
147    where Self: Sized
148    {
149        circles::angle_carve(self, radius as f64)
150    }
151
152    ///pixels on the perimeter of a cirlce with specified radius
153    fn perimimeter_pixels(&self, radius: usize) -> Vec<&Self::Data> 
154    where Self: Sized
155    {
156        let data = self.data();
157        self.perimimeter_indicies(radius)
158            .into_iter()
159            .map(|i| &data[i])
160            .collect()        
161    }
162
163    fn apply_index_mask(&self, indicies: &[usize]) -> Vec<Self::Data>
164    where Self::Data: Copy
165    {
166        indicies.iter().map(|i| self.data()[*i]).collect()
167    }
168
169    fn apply_mask<I: Image<Data=bool>>(&self, mask_image: I) -> Vec<Self::Data>
170    where Self::Data: Copy
171    {
172        let mut iter = mask_image.data().iter();
173        self.data().iter().filter(|_| *iter.next().unwrap()).map(|v| *v).collect()
174    }
175
176    fn mul<I: Image<Data=Self::Data>>(&self, rhs: I) -> OwnedImage<Self::Data>
177    where Self::Data : Mul<Self::Data, Output = Self::Data> + Copy,
178    {
179        let data = self.data().iter().zip(rhs.data().iter()).map(|(a, b)| (*a) * (*b)).collect();
180        OwnedImage::new(self.shape(), data)
181    }
182
183    fn div<I: Image<Data=Self::Data>>(&self, rhs: I) -> OwnedImage<Self::Data>
184    where Self::Data : Div<Self::Data, Output = Self::Data> + Copy,
185    {
186        let data = self.data().iter().zip(rhs.data().iter()).map(|(a, b)| (*a) / (*b)).collect();
187        OwnedImage::new(self.shape(), data)
188    }
189
190    fn map<T, F: FnMut(&Self::Data) -> T>(&self, f: F) -> OwnedImage<T>
191    {
192        let mut image : OwnedImage<_> = self.data().iter().map(f).collect();
193        image.reshape(self.shape());
194        image
195    }
196
197    #[cfg(feature = "plot")]
198    fn plot(&self, title: Option<&str>) -> ()
199    where Self::Data : DataType + Copy
200    {
201        let mut figure = Figure::new();
202        figure.axes2d()
203              .set_title(title.unwrap_or(""), &[])
204              .image(self.data().iter().map(|x| *x), self.n_rows(), self.n_cols(), None, &[]);
205        figure.show_and_keep_running().unwrap();
206    }
207
208    fn write_tiff<P: AsRef<Path>>(&self, path: P) -> Result<(), String>
209    where Self::Data : tiff_wrap::writer::ToColourType, [Self::Data]: tiff_wrap::writer::TiffValue
210    {        
211        let (width, height) = self.get_tiff_writer_shape()?;
212        let mut writer = tiff_wrap::writer::StandardTiffWriter::from_disk(path, width, height).map_err(|e| e.to_string())?;
213        writer.write_image_for(self.data()).map_err(|e| e.to_string())
214    }
215
216    fn write_big_tiff<P: AsRef<Path>>(&self, path: P) -> Result<(), String>
217    where Self::Data : tiff_wrap::writer::ToColourType, [Self::Data]: tiff_wrap::writer::TiffValue
218    {
219        let (width, height) = self.get_tiff_writer_shape()?;
220        let mut writer = tiff_wrap::writer::BigTiffWriter::from_disk(path, width, height).map_err(|e| e.to_string())?;
221        writer.write_image_for(self.data()).map_err(|e| e.to_string())
222    }
223
224    fn get_tiff_writer_shape(&self) -> Result<(u32, u32), String>
225    {
226        let (rows, cols) = self.shape();
227        let width = cols.try_into().map_err(|e: std::num::TryFromIntError| e.to_string())?;
228        let height = rows.try_into().map_err(|e: std::num::TryFromIntError| e.to_string())?;
229        Ok((width, height))
230    }
231}
232
233impl<I : Image> Image for &I
234{
235    type Data = I::Data;
236
237    fn n_rows(&self) -> usize
238    {
239        (*self).n_rows()
240    }
241
242    fn n_cols(&self) -> usize
243    {
244        (*self).n_cols()
245    }
246
247    fn data(&self) -> &[I::Data]
248    {
249        (*self).data()
250    }
251}
252
253impl<I : Image> Image for &mut I
254{
255    type Data = I::Data;
256
257    fn n_rows(&self) -> usize
258    {
259        Image::n_rows(*self)
260    }
261
262    fn n_cols(&self) -> usize
263    {
264        Image::n_cols(*self)
265    }
266
267    fn data(&self) -> &[I::Data]
268    {
269        Image::data(*self)
270    }
271}
272
273pub trait ImageMut: Image
274{
275    fn data_mut(&mut self) -> &mut [Self::Data];
276
277    /// get data with position
278    /// unsafe
279    fn get_at_mut(&mut self, row: usize, col: usize) -> &mut Self::Data
280    {
281        // println!("{row},{col}");
282        let index = self.get_index(row, col);
283        self.get_mut(index)
284    }
285
286    /// get data with index
287    /// unsafe
288    fn get_mut(&mut self, idx: usize) -> &mut Self::Data
289    {
290        &mut self.data_mut()[idx]
291    }
292
293    fn add_to<'a, I: Image<Data=Self::Data> + 'a>(&'a mut self, other: I) -> () 
294    where <Self as Image>::Data: Copy, 
295          <Self as Image>::Data: AddAssign 
296    {
297        for (value, other_value) in self.data_mut().iter_mut().zip(other.data().iter().copied())
298        {
299            *value += other_value;
300        }
301    }
302}
303
304impl<I: ImageMut> ImageMut for &mut I
305{
306    fn data_mut(&mut self) -> &mut [<Self as Image>::Data]
307    {
308        <I as ImageMut>::data_mut(*self)
309    }
310}
311
312pub trait Crop : Image
313{
314    fn take_rows(&mut self, n_rows: usize) -> Result<(), ()>;
315    fn drop_rows(&mut self, n_rows: usize) -> Result<(), ()>;
316    fn take_cols(&mut self, n: usize) -> Result<(), ()>;
317    fn drop_cols(&mut self, n_cols: usize) -> Result<(), ()>;
318
319    fn central_row_crop(&mut self, n_rows: usize) -> Result<(), ()>
320    {
321        let is_even = n_rows % 2 == 0;
322        let half = n_rows / 2;
323        let n_take = if is_even {half} else{half + 1};
324        self.take_rows(n_take).and_then(|_| self.drop_rows(half))
325    }
326
327    fn central_col_crop(&mut self, n_cols: usize) -> Result<(), ()>
328    {
329        let is_even = n_cols % 2 == 0;
330        let half = n_cols / 2;
331        let n_take = if is_even {half} else{half + 1};
332        self.take_cols(n_take).and_then(|_| self.drop_cols(half))
333    }
334
335    fn central_crop(&mut self, n_rows: usize, n_cols: usize) -> Result<(), ()>
336    {
337        self.central_row_crop(n_rows).and_then(|_| self.central_col_crop(n_cols))
338    }
339}
340
341// undefined behaviour as to what happens if you take too many rows
342// up to programmer to make sure this doesn't happen
343impl<C : Crop> Crop for &mut C
344{
345    fn take_rows(&mut self, n_rows: usize) -> Result<(), ()>
346    {
347        Crop::take_rows(*self, n_rows)
348    }
349
350    fn drop_rows(&mut self, n: usize) -> Result<(), ()>
351    {
352        Crop::drop_rows(*self, n)
353    }
354
355    fn take_cols(&mut self, n: usize) -> Result<(), ()>
356    {
357        Crop::take_cols(*self, n)
358    }
359
360    fn drop_cols(&mut self, n: usize) -> Result<(), ()>
361    {
362        Crop::drop_cols(*self, n)
363    }
364
365    fn central_row_crop(&mut self, n_rows: usize) -> Result<(), ()>
366    {
367        Crop::central_row_crop(*self, n_rows)
368    }
369
370    fn central_col_crop(&mut self, n_cols: usize) -> Result<(), ()>
371    {
372        Crop::central_col_crop(*self, n_cols)
373    }
374
375    fn central_crop(&mut self, n_rows: usize, n_cols: usize) -> Result<(), ()>
376    {
377        Crop::central_crop(*self, n_rows, n_cols)
378    }
379}
380
381pub trait Pad : Image
382{
383    fn pad_rows_start(&mut self, n_rows: usize) -> Result<(), ()>;
384    fn pad_rows_end(&mut self, n_rows: usize) -> Result<(), ()>;
385    fn pad_cols_start(&mut self, n_cols: usize) -> Result<(), ()>;
386    fn pad_cols_end(&mut self, n_cols: usize) -> Result<(), ()>;
387
388    fn central_row_pad(&mut self, n_rows: usize) -> Result<(), ()>
389    {
390        let is_even = n_rows % 2 == 0;
391        let half = n_rows / 2;
392        let start = if is_even {half} else{half + 1};
393        self.pad_rows_start(start).and_then(|_| self.pad_rows_end(half))
394    }
395
396    fn central_col_pad(&mut self, n_cols: usize) -> Result<(), ()>
397    {
398        let is_even = n_cols % 2 == 0;
399        let half = n_cols / 2;
400        let start = if is_even {half} else{half + 1};
401        self.pad_cols_start(start).and_then(|_| self.pad_cols_end(half))
402    }
403
404    fn central_pad(&mut self, n_rows: usize, n_cols: usize) -> Result<(), ()>
405    {
406        self.central_row_pad(n_rows).and_then(|_| self.central_col_pad(n_cols))
407    }
408
409    fn pad_to(&mut self, rows: usize, cols: usize) -> Result<(), ()>
410    {
411        let n_rows = rows - self.n_rows();
412        let n_cols = cols - self.n_cols();
413        self.central_pad(n_rows, n_cols)
414    }
415}
416
417impl<P : Pad> Pad for &mut P
418{
419    fn pad_rows_start(&mut self, n_rows: usize) -> Result<(), ()>
420    {
421        Pad::pad_rows_start(*self, n_rows)
422    }
423
424    fn pad_rows_end(&mut self, n: usize) -> Result<(), ()>
425    {
426        Pad::pad_rows_end(*self, n)
427    }
428
429    fn pad_cols_start(&mut self, n: usize) -> Result<(), ()>
430    {
431        Pad::pad_cols_start(*self, n)
432    }
433
434    fn pad_cols_end(&mut self, n: usize) -> Result<(), ()>
435    {
436        Pad::pad_cols_end(*self, n)
437    }
438
439    fn central_row_pad(&mut self, n_rows: usize) -> Result<(), ()>
440    {
441        Pad::central_row_pad(*self, n_rows)
442    }
443
444    fn central_col_pad(&mut self, n_cols: usize) -> Result<(), ()>
445    {
446        Pad::central_col_pad(*self, n_cols)
447    }
448
449    fn central_pad(&mut self, n_rows: usize, n_cols: usize) -> Result<(), ()>
450    {
451        Pad::central_pad(*self, n_rows, n_cols)
452    }
453
454    fn pad_to(&mut self, rows: usize, cols: usize) -> Result<(), ()>
455    {
456        Pad::pad_to(*self, rows, cols)
457    }
458}
459
460pub trait Fft
461{
462    type Output;
463    fn fft(&self) -> Self::Output;    
464    fn shifted_fft(&self) -> Self::Output;
465}
466
467impl<I> Fft for I
468where I: Image<Data=f64>
469{
470    type Output = OwnedImage<c64>;
471    fn fft(&self) -> Self::Output
472    {
473        let n_rows = self.n_rows();
474        let n_cols = self.n_cols();
475        let output_data = fft_2d(self.data(), n_rows, n_cols); //.all_data(n_rows, n_cols).into_iter().collect();
476        Self::Output::new((n_rows, n_cols), output_data)
477    }
478
479    fn shifted_fft(&self) -> Self::Output
480    {
481        let mut new_data = vec![c64::ZERO; self.size()];
482        fft_shift(self.fft().data(), &mut new_data, self.n_rows(), self.n_cols());
483        OwnedImage::new((self.n_rows(), self.n_cols()), new_data)
484    }
485}
486
487pub fn read_tiff_image_64<P: AsRef<Path>>(tiff_file: P) -> Result<OwnedImage<f64>, String>
488{
489    let mut reader = Tiff::read(tiff_file).map_err(|e| e.to_string())?;
490    let (width, height) = reader.dimensions().map_err(|e| e.to_string())?;
491    let image = reader.current_image().map_err(|e| e.to_string())?;
492    let data = image.to_f64().map_err(|e| e.to_string())?;
493    Ok(OwnedImage::new((height as usize, width as usize), data))
494}
495
496#[cfg(test)]
497mod tests 
498{
499    use super::*;
500
501    use crate::images::BorrowedImage;
502
503    fn generate(n_rows: usize, n_cols: usize) -> ((usize, usize),Vec<usize>)
504    {
505        ((n_rows, n_cols), (0..(n_rows * n_cols)).collect())
506    }
507
508    #[test]
509    fn centre_index_test() 
510    {
511        test_centre_index(4, 4, 10);
512        test_centre_index(3, 4, 6);
513        test_centre_index(3, 3, 4);
514        test_centre_index(5, 5, 12);
515    }
516
517    fn test_centre_index(n_rows: usize, n_cols: usize, expected: usize) -> ()
518    {
519        let (shape, data) = generate(n_rows, n_cols);
520        let image = BorrowedImage::new(shape, &data);
521        assert_eq!(image.centre_index(), expected);
522    }
523
524    #[test]
525    fn square_centre_test() 
526    {        
527        test_square_centre(4, 4, true);
528        test_square_centre(5, 5, false);
529    }
530
531    fn test_square_centre(n_rows: usize, n_cols: usize, expected: bool) -> ()
532    {
533        let (shape, data) = generate(n_rows, n_cols);
534        let image = BorrowedImage::new(shape, &data);
535        assert_eq!(image.square_centre(), expected);
536    }
537
538    #[test]
539    fn basic_fft()
540    {
541        let shape = (1, 4);
542        let data = vec![0.0, 1.0, 2.0, 3.0];
543        let image = BorrowedImage::new(shape, &data);
544        let fft_data = image.fft();
545        
546        let expected = [ [6.0, 0.0], [-2.0, 2.0], [-2.0, 0.0], [-2.0, -2.0] ].map(|a| c64::new(a[0], a[1]));
547        assert_eq!(fft_data.data(), expected);
548    }
549
550    #[test]
551    fn shifted_fft()
552    {
553        let shape = (1, 4);
554        let data = vec![0.0, 1.0, 2.0, 3.0];
555        let image = BorrowedImage::new(shape, &data);
556        let fft_data = image.shifted_fft();
557        
558        let expected = [ [-2.0, 0.0], [-2.0, -2.0], [6.0, 0.0], [-2.0, 2.0]].map(|a| c64::new(a[0], a[1]));
559        assert_eq!(fft_data.data(), expected);
560    }
561
562    #[test]
563    fn multiply()
564    {
565        let shape = (1, 4);
566        let data = vec![0.0, 1.0, 2.0, 3.0];
567        let image_1 = BorrowedImage::new(shape, &data);
568        let image_2 = BorrowedImage::new(shape, &data);
569        let m_image = image_1.mul(image_2);
570        let expected = [ 0.0, 1.0, 4.0, 9.0];
571        assert_eq!(m_image.data(), expected);
572    }
573
574    // #[test]
575    // fn perimenter()
576    // {
577    //     let radius = 7.0;        
578    //     let centre = (7.0, 7.0);
579    //     assert_eq!(perimeter_location(radius, centre, (0.0, 2.0)), false);
580    //     for col in 3..10
581    //     {
582    //         assert_eq!(perimeter_location(radius, centre, (0.0, col as f64)), true);
583    //     }
584    //     assert_eq!(perimeter_location(radius, centre, (0.0, 11.0)), false);
585        
586    //     assert_eq!(perimeter_location(radius, centre, (1.0, 1.0)), false);
587    //     assert_eq!(perimeter_location(radius, centre, (1.0, 2.0)), true);
588    //     assert_eq!(perimeter_location(radius, centre, (1.0, 3.0)), true);
589    //     for col in 4..9
590    //     {
591    //         assert_eq!(perimeter_location(radius, centre, (1.0, col as f64)), false);
592    //     }
593    //     assert_eq!(perimeter_location(radius, centre, (1.0, 10.0)), true);
594    //     assert_eq!(perimeter_location(radius, centre, (1.0, 11.0)), true);
595    //     assert_eq!(perimeter_location(radius, centre, (1.0, 12.0)), false);
596    // }
597
598    #[test]
599    fn perimenter_pixels()
600    {
601        let (shape, data) = generate(14, 14);
602        let image = BorrowedImage::new(shape, &data);
603        //local crossing expected
604        //let expected = [60, 61, 62, 63, 64, 65, 74, 79, 88, 93, 102, 107, 116, 121, 130, 131, 132, 133, 134, 135];
605        let expected = [61, 62, 63, 64, 65, 74, 75, 79, 80, 88, 94, 102, 108, 116, 122, 130, 131, 135, 136, 145, 146, 147, 148, 149];
606        let mut result = image.perimimeter_pixels(3);
607        result.sort();
608
609        assert_eq!(result.into_iter().map(|x| *x).collect::<Vec<usize>>().as_slice(), expected);
610    }
611
612    #[test]
613    fn drop_rows_test()
614    {
615        let (shape, data) = generate(5, 5);
616        let mut image = OwnedImage::new(shape, data);
617        assert_eq!(image.drop_rows(2).is_ok(), true);
618        assert_eq!(image.shape(), (3, 5));
619        assert_eq!(image.data(), [0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14])
620    }
621
622    #[test]
623    fn drop_cols_test()
624    {
625        let (shape, data) = generate(5, 5);
626        let mut image = OwnedImage::new(shape, data);
627        assert_eq!(image.drop_cols(2).is_ok(), true);
628        assert_eq!(image.shape(), (5, 3));
629        assert_eq!(image.data(), [0, 1, 2, 5, 6, 7, 10, 11, 12, 15, 16, 17, 20, 21, 22])
630    }
631
632    #[test]
633    fn mut_image_is_image()
634    {
635        let (shape, data) = generate(5, 5);
636        let mut image = OwnedImage::new(shape, data);
637        let image_ref = &mut image;
638        assert_eq!(image_ref.n_rows(), 5);
639        assert_eq!(image_ref.n_cols(), 5);
640        assert_eq!(image_ref.shape(), (5, 5));
641    }
642
643    #[test]
644    fn mapped_image_is_correct_shape()
645    {
646        let (shape, data) = generate(5, 5);
647        let image = OwnedImage::new(shape, data);
648        let new_image = image.map(|_| 0);
649        assert_eq!(new_image.data(), vec![0; 25]);
650        assert_eq!(new_image.shape(), image.shape())
651    }
652}
smlm_imp/lib.rs

smlm_imp/
lib.rs
