diff --git a/geomagio/ImagCDFFactory.py b/geomagio/ImagCDFFactory.py index 5d104acde4b09217fe34d3deca85731aaa22c222..3189d15a803fee175f7e2b7b0ffc7db102472fe8 100644 --- a/geomagio/ImagCDFFactory.py +++ b/geomagio/ImagCDFFactory.py @@ -15,7 +15,7 @@ from io import BytesIO import os import sys from typing import List, Optional, Union - +from datetime import datetime, timezone import numpy as np from obspy import Stream, Trace, UTCDateTime @@ -121,74 +121,80 @@ class ImagCDFFactory(TimeseriesFactory): Note: Parsing from strings is not implemented in this factory. """ raise NotImplementedError('"parse_string" not implemented') - def write_file(self, fh, timeseries: Stream, channels: List[str]): - """Write the timeseries data to a file handle in ImagCDF format. - - Parameters: - - fh: File handle to write the data. - - timeseries: ObsPy Stream containing the geomagnetic data. - - channels: List of channels to include in the output file. - """ # Create a temporary file to write the CDF data with tempfile.NamedTemporaryFile(delete=False) as tmp_file: tmp_file_path = tmp_file.name + ".cdf" try: # Initialize the CDF writer - cdf_writer = cdflib.cdfwrite.CDF(tmp_file_path, cdf_spec=None) + cdf_writer = cdflib.cdfwrite.CDF(tmp_file_path) - # Write global attributes (metadata that applies to the entire file) + # Write global attributes global_attrs = self._create_global_attributes(timeseries, channels) cdf_writer.write_globalattrs(global_attrs) - # Write time variables for each channel + # Time variables time_vars = self._create_time_stamp_variables(timeseries) for ts_name, ts_data in time_vars.items(): + # Define time variable specification var_spec = { - "Variable": ts_name, - "Data_Type": 33, # CDF TT2000 data type - "Num_Elements": 1, - "Rec_Vary": True, - "Dim_Sizes": [], - "Var_Type": "zVariable", + 'Variable': ts_name, + 'Data_Type': 33, # CDF_TIME_TT2000 + 'Num_Elements': 1, + 'Rec_Vary': True, + 'Var_Type': 'zVariable', + 'Dim_Sizes': [], + 'Sparse': 'no_sparse', + 'Compress': 6, + 'Pad': None, } - print(f"Writing time variable {ts_name} with data length: {len(ts_data)}") - cdf_writer.write_var( - var_spec=var_spec, - var_attrs=self._create_time_var_attrs(ts_name), - var_data=ts_data, - ) + # Define time variable attributes + var_attrs = self._create_time_var_attrs(ts_name) + + # Write time variable + cdf_writer.write_var(var_spec, var_attrs, ts_data) - # Write geomagnetic data variables + # Data variables for trace in timeseries: channel = trace.stats.channel var_name = f"GeomagneticField{channel}" + data_type = self._get_cdf_data_type(trace) + num_elements = 1 + CDF_CHAR = 51 + CDF_UCHAR = 52 + if data_type in [CDF_CHAR, CDF_UCHAR]: # Handle string types + num_elements = len(trace.data[0]) if len(trace.data) > 0 else 1 + var_spec = { - "Variable": var_name, - "Data_Type": self._get_cdf_data_type(trace), - "Num_Elements": 1, - "Rec_Vary": True, - "Dim_Sizes": [], - "Var_Type": "zVariable", + 'Variable': var_name, + 'Data_Type': data_type, + 'Num_Elements': num_elements, + 'Rec_Vary': True, + 'Var_Type': 'zVariable', + 'Dim_Sizes': [], + 'Sparse': 'no_sparse', + 'Compress': 6, + 'Pad': None, } - print(f"Writing data variable {var_name} with data shape: {trace.data.shape}") - cdf_writer.write_var( - var_spec=var_spec, - var_attrs=self._create_var_attrs(trace), - var_data=trace.data, - ) + + var_attrs = self._create_var_attrs(trace) + + # Write data variable + cdf_writer.write_var(var_spec, var_attrs, trace.data) + + # Close the CDF writer + cdf_writer.close() # Copy the temporary CDF file to the final file handle with open(tmp_file_path, "rb") as tmp: cdf_data = tmp.read() fh.write(cdf_data) - cdf_writer.close() - finally: # Cleanup the temporary file - print(tmp_file_path) + os.remove(tmp_file_path) + def put_timeseries( self, @@ -336,11 +342,11 @@ class ImagCDFFactory(TimeseriesFactory): 'IagaCode': {0: station}, 'ElementsRecorded': {0: ''.join(channels)}, 'PublicationLevel': {0: publication_level}, - 'PublicationDate': {0: UTCDateTime.now().strftime("%Y-%m-%dT%H:%M:%SZ")}, + 'PublicationDate': {0: [cdflib.cdfepoch.timestamp_to_tt2000(datetime.timestamp(datetime.now(timezone.utc))), "cdf_time_tt2000"]}, 'ObservatoryName': {0: observatory_name}, - 'Latitude': {0: latitude}, - 'Longitude': {0: longitude}, - 'Elevation': {0: elevation}, + 'Latitude': {0: np.array([latitude], dtype=np.float64)}, + 'Longitude': {0: np.array([longitude], dtype=np.float64)}, + 'Elevation': {0: np.array([elevation], dtype=np.float64)}, 'Institution': {0: institution}, 'VectorSensOrient': {0: vector_orientation}, #remove F - because its a calculation, not an element? 'StandardLevel': {0: 'None'}, # Set to 'None' @@ -436,7 +442,7 @@ class ImagCDFFactory(TimeseriesFactory): return var_spec def _create_var_attrs(self, trace: Trace) -> dict: - print(trace.stats) + # print(trace.stats) channel = trace.stats.channel fieldnam = f"Geomagnetic Field Element {channel}" # “Geomagnetic Field Element †+ the element code or “Temperature †+ the name of the location where the temperature was recorded. units = '' # Must be one of “nTâ€, “Degrees of arc†or “Celsiusâ€