Madrigal HDF5 support, memory-safe loader, and a working GOES data source - #3
Open
w2naf wants to merge 11 commits into
Open
Madrigal HDF5 support, memory-safe loader, and a working GOES data source#3w2naf wants to merge 11 commits into
w2naf wants to merge 11 commits into
Conversation
Replaces the prototype feature_hdf5 loader with a working one: - New gen_lib.load_madrigal_hdf5() using pytables (h5py 3.15 cannot decode Madrigal's float precision spec). - Per-file ut1_unix offset detection + correction; older ham files carry a local-timezone offset (empirically: most of 2023 and 2024 are +5h/+6h; other years clean). - Widened band table to cover 2200m LF, WARC bands, 2m VHF; out-of- band spots now flagged with -1 instead of 0. - ssrc (PSK/RBN/WSP) mapped to harc_plot's existing integer source codes (1=WSPRNet, 2=RBN, 3=PSKReporter). - Uses Madrigal-provided latcen/loncen for md_lat/md_long instead of recomputing midpoints. - Fixed regional_filter(kind='endpoints'): restored logical_or (prior pd.concat implementation duplicated rows where both endpoints were inside the region). - Re-enabled regional_filter call in load_spots_csv (was commented out). - Removed dask dependency (no real speedup benefit; complicated memory profile on large files). - Added guard around loc_sources filter for HDF5 path (CSV-only fields). - Deleted calculate_histograms_h5.py: redundant with calculate_histograms.py now that load_spots_csv transparently handles HDF5. - Cleanup: removed scratch notebook nc_test.ipynb; restored sample 20170701.data.nc.bz2; added tables to requirements.txt. Validated end-to-end on 2017-07-01 (clean, 3.2M rows) and 2023-06-21 (offset-corrected, 114M rows); full histogram pipeline produces a valid 125 KB NetCDF at 20s/day on a clean file. Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
…value convention - Replace 2017 CSV sample (22 MB, pre-PSKReporter) with data/madrigal/rsd2015-06-10.01.hdf5 (39 MB, all three networks, ~963k rows). First day with RBN+WSPRNet+PSKReporter coverage that stays under GitHub's 50 MB soft limit. - Update scripts/basic_histogram/histograms.py to the new date and set data_sources=[1,2,3] explicitly so PSK is not filtered out. - Fix band value convention in load_madrigal_hdf5: the prior commit emitted int(MHz) keys (1,3,7,14,21,28), but calculate_histograms filters via df["band"] == band["meters"] expecting wavelength values. This silently mis-bucketed all 30m WSPR traffic as the only 10m band result (band=10 int(MHz) happened to coincide with the 10m-meters filter value). Change _BAND_EDGES_MHZ to emit meters values to match the historical convention. - Verified end-to-end: 753k spots across all 6 bands, PSK 45% / RBN 18% / WSP 37%, correct diurnal patterns, per-band sidebar maps now populate correctly (the N=0 issue in the prior plots was a symptom of the same bug, not a separate issue). - Regenerate data/histograms/World/20150610.data.nc.bz2 (635 KB) and example_plots/20150610.*.png. - Rewrite README for Python 3.10+, conda-based install (cartopy + pytables via conda-forge), new data-input section, ut1_unix local-timezone offset note. Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
Prep for the 5.1 TB solar-cycle plotting pass: fix the memory envelope so multi-hundred-million-row Madrigal days process without OOM. Root cause: the 78 GB RSS previously observed on 114M-row days came from string-column overhead (call_sign_tx, call_sign_rx, smode loaded as Python objects) pulled in by the full recarray read. Numeric data for the same file is ~3 GB when packed into compact dtypes. Changes: - load_madrigal_hdf5: per-field pytables reads (no full recarray materialization). Compact dtypes: float32 for coords / distances, int8 for source, int16 for band. New `columns` kwarg for selective loading with an explicit whitelist + error on unknown names. - load_spots_csv: threads a minimal column set into the HDF5 loader based on filter_region_kind. mids-only (the default) skips tx/rx coords. - Vectorize ut_hrs / slt_mid via the pandas .dt accessor (was a per-row Python lambda — multi-minute cost on 100M-row days). - Replace df.source.map(lambda x: x in ...) with .isin(...). Benchmark on rsd2023-06-21.01.hdf5 (114M rows): Before: 78 GB RSS, process killed After: 2.89 GB RSS, 373s wall, 2.63 GB DataFrame Regression-tested basic_histogram on 2015-06-10: per-band spot counts identical to pre-refactor output. The 20150610.data.nc.bz2 fixture is regenerated (identical values, new dtypes in the underlying columns). Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
scripts/solarcycle/batch_histograms.py runs harc_plot's per-day histogram pipeline across a date range using multiprocessing.Pool, producing one YYYYMMDD.data.nc.bz2 per day under a config-keyed output dir. Key design choices: - spawn context avoids fork-safety issues with pytables open HDF5 handles. - The driver clears the output dir once up-front (honoring --reprocess), then workers always run reprocess=False. This fixes a race where per-worker reprocess=True was clearing the output dir concurrently and non-deterministically truncating results. - --network selects WSPRNet / RBN / PSKReporter / All, for the three-networks-separate workflow. - imap_unordered + chunksize=1 + tqdm gives a responsive progress bar that tracks actual finished days (load balancing across highly variable day sizes — 2009 files are sub-MB, 2023+ files are 100+ MB). Validation: - Serial-vs-parallel byte-identical on 2015-06-10..14 (5 days). - 4-worker parallel on 5-day range: 18s wall vs 52s serial (~3x). - End-to-end on the 114M-row 2023-06-21 day: 7:03 wall, 11.4 GB peak RSS, valid 743 KB NetCDF with correct per-band distribution (20m 40.8M, 15m 20.5M, 40m 8.9M — plausible June 2023). - Extrapolated full-archive budget at 32 workers on a 32-core / 566 GB host: ~8-22 hrs per network depending on year mix. Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
Last night's WSPRNet run crashed after 5 hrs / 3040 successful outputs when a worker hit a truncated 2023 Madrigal HDF5 file: the unhandled exception propagated through imap_unordered and killed the whole pool. One bad file shouldn't abort the whole pipeline. Changes: - Wrap _process_day body in try/except; workers now return (date, None | error_message) tuples. - Driver collects failures, prints a summary, and writes _failed_dates.txt alongside the outputs for later triage. A file scan turned up 21 bad files out of 6118 (0.34%): - 19x HDF5ExtError (truncated / corrupt). - 2x NoSuchNodeError (missing /Data/Table Layout). - Notable cluster: 14 consecutive days 2025-02-07..20, likely an upstream ingestion outage at MIT Haystack. These are upstream data issues, not pipeline bugs. The hardened script skips them and continues. Relaunched WSPRNet at 10:59 UTC; it resumes by skipping the 3040 existing outputs via the reprocess=False existence check. Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
scripts/solarcycle/plot_solarcycle.py aggregates the per-day NetCDFs produced by batch_histograms.py into a single multi-year panel plot, one per network. Adapted from the 2020 agu2020_plot_* scripts with updates for the new output layout, an argparse interface, and per-network start-date defaults that reflect when each network enters Madrigal (WSPRNet/PSKReporter: 2015-01-01; RBN: 2009-02-25). Concatenated xarray Dataset is pickled to a cache so iterating on plot styling does not re-read the ~4000 daily NetCDFs each time. Validated end-to-end: - 3-month WSPRNet slice (N=55.8M, all 6 bands, correct diurnal patterns). - Full 11-year WSPRNet archive (N=7.43B across 3999 per-day NetCDFs): solar-cycle panel plot with F10.7 overlay showing Cycle 24 minimum and Cycle 25 ramp clearly. Also cleaned up an ipdb breakpoint that had been left at the end of the 2020 template. Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
- Re-ran data/omni/update_omni_data.py (unchanged script, already in repo) to pull OMNI2 hourly through 2026-03-21. The prior file stopped at 2020 day 310 / Nov 2020, missing the entire Cycle 25 ramp. New file: 229,848 rows, verified loadable via harc_plot.omni.Omni. - Added a "Data sources and acknowledgements" section to README.md with a paste-ready blockquote citing NASA OMNIWeb / OMNI2 and the DRAO Penticton / NRC Canada origin of the F10.7 measurement. - Added a module-level docstring to update_omni_data.py carrying the same attribution, plus a note that the other bundled indices (Kp, Dst, ap, AE/AL/AU) have their own attribution chain through ISGI / GFZ Potsdam / Kyoto WDC. Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
Two backwards-compatible fixes to the Kp/ring-current panel, found while validating storm-day plots against intense events. omni.py: plot_dst_kp hardcoded set_ylim(-200,100), which silently clipped intense storms off the axis entirely (-223 nT on 2015-03-17, -412 nT on 2024-05-10). Limits are now expand-only, using (-200,100) as a floor rather than a cap. They are computed over the plotted window only: the SYM-H branch hands matplotlib the full multi-decade series and relies on set_xlim to window it, so a naive min/max would pick up the all-time extreme. Plots whose data already fits within the old limits are unchanged. visualize_histograms.py: the Kp/SYM-H panel hardcoded dst_param='SYM-H', but SymH defaults to years=[2016,2017] and no Kyoto SYM file exists past 2019, so any plot outside 2016-2017 rendered an empty geomagnetic panel with no warning. dst_param is now a plot() kwarg with the default unchanged; callers spanning arbitrary epochs can select 'Dst_nT' (hourly OMNI, covers the whole archive). Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
…ble panel order The GOES X-ray panel has been silently empty for every caller. read_goes() downloads from ftp://satdat.ngdc.noaa.gov, which NOAA retired along with anonymous FTP, so it always hits its bare `except`, prints 'GOES Data ERROR.', and returns None. Goeser then skips the satellite and plots a blank axis. It also only ever targeted GOES-13/15, which never covered 2021 onward. read_goes_ncei() replaces it against the two NCEI HTTPS archives, which between them span 2009 to present: pre-R (GOES 8-15) GOES Space Environment Monitor, one netCDF per month, variables A_AVG / B_AVG, time_tag GOES-R (GOES 16-19) L2 xrsf-l2-avg1m, one netCDF per day, variables xrsa_flux / xrsb_flux, time Both are normalized to the A_AVG / B_AVG names goes_plot() and find_flares() already expect, so nothing downstream changes. Downloads go to a .part file and are renamed, so an interrupted fetch cannot leave a truncated file that a later run treats as a valid cache hit. Corrupt cached files are dropped and re-fetched. Flare magnitudes are NOT comparable across satellite generations. For the 2017-09-06 flare: GOES-15 peaks at 9.33e-4 W/m^2 (X9.3, the accepted classification), GOES-13 at 1.04e-3 (X10.4), GOES-16 -- then still in post-launch test -- at 5.17e-4 (X5.2). Nearly 2x for one event. goes_xrs_candidates() therefore prefers the primary operational satellite for each epoch, which reproduces the standard classifications. Also fixes a pre-existing crash in find_flares(): `idx_max is np.nan` never matches pd.NaT, so an all-NaN window (eclipse, data gap, flagged-bad stretch) put NaT into the key list and the b_avg[keys] lookup raised "KeyError: '[NaT] not in index'". Windows are now tested for emptiness up front, which also avoids the deprecated all-NA idxmax call. Goeser gains a 120-day span guard -- 1-minute flux is meaningless on a multi-year axis and fetching it would pull thousands of daily files -- so the solar-cycle plots now skip the panel deliberately instead of by accident. An empty panel is annotated "No GOES X-ray data available for this interval" so it is not misread as a flare-free period. Finally, the three context panels (Kp/SYM-H, GOES, F10.7) become closures driven by a new context_panels parameter, so callers choose their order and the panel letters follow. The default is unchanged, keeping existing figures lettered as before. Invalid or duplicate keys raise before the no-data early return, so a typo fails on every call rather than only on calls with data. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
2 participants
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
Brings
feature_hdf5intomaster. This is a fast-forward —masterhas no commits that this branch lacks, so there is nothing to resolve.The headline is Madrigal HDF5 support:
mastercannot read the current data format at all, and the memory-safe loader here is what makes a full-archive pass feasible.What's in it
rsdYYYY-MM-DD.01.hdf5directly; auto-detect HDF5 vs. legacy CSVut1_unixlocal-timezone offset against the authoritativey/m/d/h/m/scolumnsscripts/solarcycle/batch_histograms.py— parallel per-day processing with per-day exception isolationscripts/solarcycle/plot_solarcycle.py— multi-year aggregation with F10.7 overlaydst_paramBehavior changes worth knowing about
These are not obvious from the commit titles:
plot_goes=Truenow performs network I/O. It previously always failed —read_goes()targetsftp://satdat.ngdc.noaa.gov, which NOAA retired with anonymous FTP, so it hit its bareexceptand returnedNone, and every caller silently got a blank panel. It now fetches from NCEI and caches underdata/goes/. Bounded by a 120-day span guard, so multi-year plots skip the panel rather than pulling thousands of daily files.GOES flare magnitudes are not comparable across satellite generations. For 2017-09-06: GOES-15 peaks at 9.33e-4 W/m² (X9.3, the accepted classification), GOES-13 at 1.04e-3 (X10.4), GOES-16 at 5.17e-4 (X5.2) — nearly 2× for one event.
goes_xrs_candidates()prefers the primary operational satellite per epoch, which reproduces standard classifications. Anyone comparing X-class numbers across the archive needs to account for this.Band value convention fix (
b8b0bad). This changes data semantics, so histogram NetCDFs generated before it may not be comparable to ones generated after.Ring-current panel defaults are unchanged, but note
dst_param='SYM-H'only has Kyoto data for the yearsSymHis configured with, and no Kyoto file exists past 2019. Callers spanning arbitrary epochs should pass'Dst_nT'(hourly OMNI, full archive).context_panelsmakes the Kp/GOES/F10.7 panel order caller-selectable. The default is unchanged, so existing figures keep their panel lettering.Also fixed
find_flares()crashed on any all-NaN window:idx_max is np.nannever matchespd.NaT, so NaT entered the key list andb_avg[keys]raisedKeyError: '[NaT] not in index'..ipynb_checkpointsfiles.Things to weigh
rsd2015-06-10.01.hdf5, now covering all three networks), plus a refreshed OMNI bundle and regenerated example PNGs.wave_search_naf(63 commits ahead of master),wave_search(18), andagu2020_solarcycle(18) are all based onmasterand will each pick up 10 commits. Worth coordinating if any are active.Validation
Sample workflow (
scripts/basic_histogram/histograms.py) runs clean on the bundled Madrigal sample. The full-archive pass has been exercised on the 11-year WSPRNet record (N = 7.43e9 spots, 3999 daily NetCDFs). GOES verified across both archive eras — GOES-15 for 2015 and 2017, GOES-18 for 2024 — with flare classifications matching accepted values.🤖 Generated with Claude Code