Merge pull request #86 from MWATelescope/pyo3_stub_chrono_fix

Pyo3 stub chrono fix-Reentract check-BSCALE

CHANGELOG.md

@@ -2,6 +2,16 @@
 
 Changes in each release are listed below.
 
+## 1.8.0 11-Nov-2024
+
+* mwalib now will detect and raise an error (`MwalibError::Fits.CfitsioIsNotReentrant`) if the CFITSIO library that mwalib is linked with has been built without the `-D_REENTRANT` directive (github issue #82).
+* Expose the FITS BSCALE in image HDUs as a single value `bscale` in `CorrelatorContext` (github issue #85).
+  * For most Legacy MWA observations, this may be a value other than 1.0. For MWAX correlator this will always be 1.0.
+  * This is mainly of interest to EoR researchers who are trying to implement Van Vleck corrections.
+* Python type stubs:
+  * Fixed pyo3 decorators to allow stub_gen to work properly / generate python stubs correctly.
+  * Used prerelease version of pyo3_stub_gen to ensure Chrono::FixedTimeOffset can have a stub generated in `MetafitsContext`.
+
 ## 1.7.2 8-Nov-2024
 
 * Update to fitsio 0.21.6 - provides FitsFile.file_path() method.
@@ -14,7 +24,7 @@ Changes in each release are listed below.
 ## 1.7.0 23-Oct-2024 (Yanked- MWAFitsFile is no longer to be used)
 
 * Bumped MSRV to 1.65.
-* Update fitsio to 0.21 and fitsio-sys to 0.5.  To make v0.21 work, a new struct, MWAFitsFile is used inplace of FitsFile, as FitsFile no longer carries the `filename` property which is needed by mwalib.
+* Update fitsio to 0.21 and fitsio-sys to 0.5. To make v0.21 work, a new struct, MWAFitsFile is used inplace of FitsFile, as FitsFile no longer carries the `filename` property which is needed by mwalib.
 * Removed Rust Report Card from README status badges. Looks like this service is abandonded.
 * Added Python .pyi stub generation to provide mwalib Python users with type and docstring information. The mwalib.pyi should get baked into the python wheels released to github and Pypi. See `bin/README.md` for caveats and more details.
 * Added CI to test compilation against cfitsio 3.x and 4.x when not using the `cfitsio-static` feature.

Cargo.toml

@@ -1,6 +1,6 @@
 [package]
 name = "mwalib"
-version = "1.7.2"
+version = "1.8.0"
 homepage = "https://github.com/MWATelescope/mwalib"
 repository = "https://github.com/MWATelescope/mwalib"
 readme = "README.md"
@@ -54,7 +54,9 @@ env_logger = { version = "~0.10", optional = true }
 ndarray = { version = "~0.16", optional = true }
 numpy = { version = "~0.22", optional = true }
 pyo3 = { version = "~0.22", features = ["chrono", "extension-module", "macros"], optional = true }
-pyo3-stub-gen = { version = "~0.6", optional = true }
+# pyo3-stub-gen = { version = "~0.6", optional = true }
+# The below pyo3-stub-gen branch has not been merged with main as of 8-Nov-2024. once it has been it should be fine to go back to using the main branch again.
+pyo3-stub-gen = { branch = "chrono_feature", optional = true, git = "https://github.com/ThatOneShortGuy/pyo3-stub-gen" }
 pyo3-stub-gen-derive = { version = "~0.6", optional = true }
 
 # "examples" feature.

bin/stub_gen.rs

@@ -39,23 +39,6 @@ fn fix_stubs() -> anyhow::Result<()> {
     // After the stub is generated, we have some "manual" fixes to do
     let stubfile = String::from("mwalib.pyi");
 
-    // Import datetime
-    insert_stub_below(&stubfile, "import typing\n", "import datetime\n")?;
-
-    // Add sched_start_utc (Chrono::DateTime<FixedOffset> is not supported yet)
-    insert_stub_below(
-        &stubfile,
-        "    sched_end_unix_time_ms: int\n",
-        "    sched_start_utc: datetime.datetime\n",
-    )?;
-
-    // Add sched_end_utc (Chrono::DateTime<FixedOffset> is not supported yet)
-    insert_stub_below(
-        &stubfile,
-        "    sched_start_utc: datetime.datetime\n",
-        "    sched_end_utc: datetime.datetime\n",
-    )?;
-
     // Replace the constructors as pyo3_stub_gen seems to ignore the text_signature
     replace_stub(&stubfile,"def __new__(cls,metafits_filename,mwa_version = ...): ...","def __new__(cls, metafits_filename: str, mwa_version: typing.Optional[MWAVersion]=None)->MetafitsContext:\n        ...\n",)?;
 
@@ -87,6 +70,8 @@ fn fix_stubs() -> anyhow::Result<()> {
 /// Inserts new items in the stubfile for when pyo3 cant do it
 /// properly.
 ///
+/// Currently this function is not used, but kept in case it is needed later!
+///
 /// This will:
 /// * Open the created stubfile
 /// * Find the line in `string_to_find` (fails if not found)
@@ -106,6 +91,7 @@ fn fix_stubs() -> anyhow::Result<()> {
 /// * Result Ok if stub file was modified successfully.
 ///
 ///
+#[allow(dead_code)]
 #[cfg(feature = "python")]
 fn insert_stub_below<P: AsRef<Path>>(
     stubfile: P,

examples/build_c_examples.sh

@@ -2,7 +2,7 @@
 
 set -eux
 
-cargo build --release
+cargo build --release --features=cfitsio-static
 
 # Compile the example C code
 gcc -O3 \

mwalib.pyi

@@ -332,6 +332,7 @@ class CorrelatorContext:
     num_gpubox_files: int
     gpubox_batches: list[GpuBoxBatch]
     gpubox_time_map: dict[int, dict[int, tuple[int, int]]]
+    bscale: float
     def __new__(cls, metafits_filename: str, gpubox_filenames: list[str]) -> CorrelatorContext: ...
     def get_fine_chan_freqs_hz_array(self, corr_coarse_chan_indices: typing.Sequence[int]) -> list[float]:
         r"""

src/antenna/mod.rs

@@ -7,15 +7,16 @@
 use crate::rfinput::*;
 use std::fmt;
 
+#[cfg(feature = "python")]
+use pyo3::prelude::*;
 #[cfg(feature = "python")]
 use pyo3_stub_gen_derive::gen_stub_pyclass;
 
 #[cfg(test)]
 mod test;
 
 /// Structure for storing MWA antennas (tiles without polarisation) information from the metafits file
-#[cfg_attr(feature = "python", gen_stub_pyclass)]
-#[cfg_attr(feature = "python", pyo3::pyclass(get_all, set_all))]
+#[cfg_attr(feature = "python", gen_stub_pyclass, pyclass(get_all, set_all))]
 #[derive(Clone)]
 pub struct Antenna {
     /// This is the antenna number.

src/baseline/mod.rs

@@ -7,14 +7,16 @@
 use crate::misc;
 use std::fmt;
 
+#[cfg(feature = "python")]
+use pyo3::prelude::*;
 #[cfg(feature = "python")]
 use pyo3_stub_gen_derive::gen_stub_pyclass;
 
 #[cfg(test)]
 mod test;
+
 /// This is a struct for our baselines, so callers know the antenna ordering
-#[cfg_attr(feature = "python", gen_stub_pyclass)]
-#[cfg_attr(feature = "python", pyo3::pyclass(get_all, set_all))]
+#[cfg_attr(feature = "python", gen_stub_pyclass, pyclass(get_all, set_all))]
 #[derive(Clone)]
 pub struct Baseline {
     /// Index in the mwalibContext.antenna array for antenna1 for this baseline    

src/coarse_channel/mod.rs

@@ -12,15 +12,16 @@ use error::CoarseChannelError;
 use fitsio::{hdu::FitsHdu, FitsFile};
 use std::fmt;
 
+#[cfg(feature = "python")]
+use pyo3::prelude::*;
 #[cfg(feature = "python")]
 use pyo3_stub_gen_derive::gen_stub_pyclass;
 
 #[cfg(test)]
 mod test;
 
 /// This is a struct for coarse channels
-#[cfg_attr(feature = "python", gen_stub_pyclass)]
-#[cfg_attr(feature = "python", pyo3::pyclass(get_all, set_all))]
+#[cfg_attr(feature = "python", gen_stub_pyclass, pyclass(get_all, set_all))]
 #[derive(Clone)]
 pub struct CoarseChannel {
     /// Correlator channel is 0 indexed (0..N-1)

src/convert/mod.rs

@@ -12,6 +12,11 @@ use crate::rfinput::*;
 use log::trace;
 use std::fmt;
 
+#[cfg(feature = "python")]
+use pyo3::prelude::*;
+#[cfg(feature = "python")]
+use pyo3_stub_gen_derive::gen_stub_pyclass;
+
 #[cfg(test)]
 mod test;
 
@@ -35,8 +40,8 @@ fn fine_pfb_reorder(input: usize) -> usize {
 }
 
 /// Structure for storing where in the input visibilities to get the specified baseline when converting
+#[cfg_attr(feature = "python", gen_stub_pyclass, pyclass(get_all, set_all))]
 #[derive(Clone)]
-#[cfg_attr(feature = "python", pyo3::pyclass(get_all, set_all))]
 pub(crate) struct LegacyConversionBaseline {
     pub baseline: usize,    // baseline index
     pub ant1: usize,        // antenna1 index

src/correlator_context/mod.rs

@@ -5,6 +5,7 @@
 //! The main interface to MWA data.
 
 use fitsio::FitsFile;
+use log::warn;
 use std::collections::BTreeMap;
 use std::fmt;
 use std::path::Path;
@@ -17,9 +18,10 @@ use crate::metafits_context::*;
 use crate::timestep::*;
 use crate::*;
 
+#[cfg(feature = "python")]
+use pyo3::prelude::*;
 #[cfg(feature = "python")]
 use pyo3_stub_gen_derive::gen_stub_pyclass;
-
 #[cfg(feature = "python")]
 mod python;
 
@@ -29,8 +31,7 @@ mod test;
 ///
 /// This represents the basic metadata and methods for an MWA correlator observation.
 ///
-#[cfg_attr(feature = "python", gen_stub_pyclass)]
-#[cfg_attr(feature = "python", pyo3::pyclass(get_all, set_all))]
+#[cfg_attr(feature = "python", gen_stub_pyclass, pyclass(get_all, set_all))]
 #[derive(Debug)]
 pub struct CorrelatorContext {
     /// Observation Metadata obtained from the metafits file
@@ -124,6 +125,8 @@ pub struct CorrelatorContext {
     pub gpubox_time_map: BTreeMap<u64, BTreeMap<usize, (usize, usize)>>,
     /// A conversion table to optimise reading of legacy MWA HDUs
     pub(crate) legacy_conversion_table: Vec<LegacyConversionBaseline>,
+    /// BSCALE- FITS BSCALE or SCALEFAC value set on the visibility HDUs (used in Legacy Correlator only)
+    pub bscale: f32,
 }
 
 impl CorrelatorContext {
@@ -147,6 +150,11 @@ impl CorrelatorContext {
         metafits_filename: P,
         gpubox_filenames: &[P2],
     ) -> Result<Self, MwalibError> {
+        // Check CFITSIO is reentrant before proceeding
+        if !fits_read::is_fitsio_reentrant() {
+            return Err(MwalibError::Fits(FitsError::CfitsioIsNotReentrant));
+        }
+
         Self::new_inner(metafits_filename.as_ref(), gpubox_filenames)
     }
 
@@ -167,6 +175,15 @@ impl CorrelatorContext {
         // Update the metafits now that we know the mwa_version
         metafits_context.mwa_version = Some(gpubox_info.mwa_version);
 
+        // Determine BSCALE from 1st visibility HDU of each gpubox file
+        let bscale: f32 = match metafits_context.mwa_version {
+            Some(MWAVersion::CorrOldLegacy) | Some(MWAVersion::CorrLegacy) => {
+                Self::get_bscale_from_gpubox_files(metafits_context.obs_id, gpubox_filenames)
+            }
+            // Only Legacy MWA Correlator observations have BSCALE set to anything other than 1.0
+            _ => 1.0,
+        };
+
         // Populate metafits coarse channels and timesteps now that we know what MWA Version we are dealing with
         // Populate the coarse channels
         metafits_context.populate_expected_coarse_channels(gpubox_info.mwa_version)?;
@@ -390,6 +407,7 @@ impl CorrelatorContext {
             num_timestep_coarse_chan_weight_floats: weight_floats,
             num_gpubox_files: gpubox_filenames.len(),
             legacy_conversion_table,
+            bscale,
         })
     }
 
@@ -902,6 +920,63 @@ impl CorrelatorContext {
 
         Ok(())
     }
+
+    /// Returns the BSCALE/SCALEFAC value used in the visibility FITS HDUs. The value
+    /// is only non-one in Legacy Correlator observations and is most of interest
+    /// only to EoR researchers.
+    ///    
+    /// # Arguments
+    ///
+    /// * `obs_id` - The observation ID of this observation
+    ///
+    /// * `gpubox_filenames` - A reference to a vector of paths to the gpubox files
+    ///
+    /// # Returns
+    ///
+    /// * an f32 which is either the BSCALE/SCALEFAC present in the visibility gpuboxfiles or 1.0
+    ///
+    fn get_bscale_from_gpubox_files<P: AsRef<Path>>(obs_id: u32, gpubox_filenames: &[P]) -> f32 {
+        let mut scale_factor: Option<f32> = None;
+
+        for raw_path in gpubox_filenames {
+            // Open each gpubox file
+            let mut fptr = fitsio::FitsFile::open(raw_path).unwrap();
+            // Only check the first visibilities HDU in each gpubox file
+            let hdu1 = fits_open_hdu!(&mut fptr, 1).unwrap();
+            // Check for SCALEFAC and BSCALE
+            for key in ["SCALEFAC", "BSCALE"] {
+                let this_scale_factor: Option<f32> =
+                    get_optional_fits_key!(&mut fptr, &hdu1, key).unwrap();
+                match (scale_factor, this_scale_factor) {
+                    (Some(sf), Some(this_sf)) => {
+                        assert!(
+                            ((sf - this_sf).abs() < f32::EPSILON),
+                            "Different scale factors found in gpubox files: {} and {}",
+                            sf,
+                            this_sf
+                        );
+                    }
+                    (None, Some(this_sf)) => {
+                        scale_factor = Some(this_sf);
+                    }
+                    _ => {}
+                }
+            }
+        }
+        // If BSCALE/SCALEFAC was found then return it, otherwise
+        // Check if the obsid is < 1096160568 then return 0.25 otherwise return 1.0
+        // TODO we should really check this "0.25" assumption!
+        match (scale_factor, obs_id) {
+            (Some(sf), _) => sf,
+            // according to pyuvdata "correlator did a divide by 4 before october 2014"
+            // https://github.com/RadioAstronomySoftwareGroup/pyuvdata/blob/05ee100af2e4e11c9d291c9eafc937578ef01763/src/pyuvdata/uvdata/mwa_corr_fits.py#L1464
+            (None, t) if t < 1096160568 => 0.25,
+            _ => {
+                warn!("No scale factor found in metafits or gpufits files, defaulting to 1.0");
+                1.0
+            }
+        }
+    }
 }
 
 /// Implements fmt::Display for CorrelatorContext struct
@@ -958,6 +1033,8 @@ impl fmt::Display for CorrelatorContext {
             Memory usage per scan:      {scan_size} MiB,
 
             gpubox batches:             {batches:#?},
+
+            BSCALE/SCALEFAC:            {bscale},
         )"#,
             metafits_context = self.metafits_context,
             corr_ver = self.mwa_version,
@@ -992,6 +1069,7 @@ impl fmt::Display for CorrelatorContext {
             hdu_size = size,
             scan_size = size * self.num_gpubox_files as f64,
             batches = self.gpubox_batches,
+            bscale = self.bscale
         )
     }
 }

src/correlator_context/python.rs

@@ -15,8 +15,6 @@ use numpy::PyArray2;
 #[cfg(feature = "python")]
 use numpy::PyArray3;
 #[cfg(feature = "python")]
-use pyo3::prelude::*;
-#[cfg(feature = "python")]
 use pyo3_stub_gen_derive::gen_stub_pymethods;
 
 #[cfg_attr(feature = "python", gen_stub_pymethods)]

src/correlator_context/test.rs

@@ -177,6 +177,8 @@ fn test_context_legacy_v1() {
         context.metafits_context.num_metafits_fine_chan_freqs,
         context.metafits_context.metafits_fine_chan_freqs_hz.len()
     );
+
+    assert_eq!(context.bscale, 0.5);
 }
 
 #[test]
@@ -272,7 +274,9 @@ fn test_context_mwax() {
     assert_eq!(
         context.num_timestep_coarse_chan_weight_floats,
         ((128 * 129) / 2) * 4
-    )
+    );
+
+    assert_eq!(context.bscale, 1.0);
 }
 
 #[test]

src/error.rs

@@ -9,6 +9,10 @@ use thiserror::Error;
 /// MwalibError subtypes
 #[derive(Error, Debug)]
 pub enum MwalibError {
+    /// An error derived from `CfitsioError`.
+    //#[error("{0}")]
+    //Cfitsio(#[from] crate::fits_read::error::CfitsioError),
+
     /// An error derived from `FitsError`.
     #[error("{0}")]
     Fits(#[from] crate::fits_read::error::FitsError),