Analyze attributed climate observations

// Actual NOAA NCEI GHCN-Daily observations, not simulated measurements.
// USW00094728 / NY CITY CNTRL PARK, all 366 calendar days of leap year 2024.
// Pinned source records, scaling rules, citations and hashes accompany the CSV.
// Run from the repository root; execution is offline and needs no map packs.
import aner.data
import aner.plot
import aner.units

// Acquisition was performed once, outside this analysis. The offline converter
// preserves flags, scales tenths to degC/mm and maps -9999 to explicit NA.
// Inference is a separate source pass; review it before accepting unfamiliar data.
let source = data.scan_csv_auto("examples/fixtures/climate_noaa_2024/observations.csv",
    missing: "NA", memory_bytes: 8388608, sample_rows: 0)
let observations = source.collect(max_rows: 0, max_bytes: 8388608, batch_rows: 64)
print("BASIC INSPECTION: cached structure and descriptive statistics")
let basic = observations.overview()
print(basic.summary())
print("Basic inspection is available without a deep profile")
print(!observations.has_profile() && basic.rows() == 366)
print(basic.missing("TMAX_DegC") == 0 && basic.missing("PRCP_mm") == 0)
print(observations.complete() && observations.origin_complete())
print(observations.column_type(2) == "String")

// Decide to inspect distributions only after seeing the basic structure.
// This deliberate deep pass reads the immutable snapshot, not the CSV again.
// If deep diagnostics are known to be needed before loading, use
// source.with_profile(512, 32).collect(...) and table.profile() instead.
let deep = observations.select(List("TMAX_DegC", "TMIN_DegC", "PRCP_mm", "PRCP_mflag"))
    .with_profile_workers(4)
    .profile(max_bytes: 8388608, batch_rows: 64, sample_rows: 512, max_categories: 32)
print("DEEP PROFILE: exact observed quantiles for this small snapshot")
print(deep.summary())
print(deep.quantiles_exact("TMAX_DegC") && deep.sample_count("TMAX_DegC") == 366)
print("Explicit native CPU profile workers")
print(deep.profile_workers())
print(deep.minimum_text("TMAX_DegC"))
print(deep.maximum_text("TMAX_DegC"))
print(deep.quantile_text("TMAX_DegC", q: 0.5))

// A blank NOAA quality flag means no listed quality test failed; it is not a
// proof of perfect measurement. Empty flags are present Strings, not NA values.
let quality_clear = (data.col("TMAX_qflag") == "") &&
    (data.col("TMIN_qflag") == "") && (data.col("PRCP_qflag") == "")
let complete_measurements = data.col("TMAX_DegC").is_present() &&
    data.col("TMIN_DegC").is_present() && data.col("PRCP_mm").is_present()
let physically_ordered = (data.col("TMAX_DegC") >= data.col("TMIN_DegC")) &&
    (data.col("PRCP_mm") >= 0.0@mm)
let accepted = observations.filter(quality_clear && complete_measurements && physically_ordered)
    .collect(max_rows: 0, max_bytes: 8388608, batch_rows: 64)
print("QC: retained rows and source immutability")
print(accepted.rows())
print(observations.rows() == 366 && accepted.rows() == 366)

// Keep trace information. T means trace, even when its numeric total is zero;
// neither replacing NA with zero nor calling every recorded zero a dry day is
// scientifically justified by this profile alone. No imputation occurs here.
let trace = accepted.filter(data.col("PRCP_mflag") == "T")
    .summarize(max_bytes: 8388608, batch_rows: 64)
let recorded_positive = accepted.filter(data.col("PRCP_mm") > 0.0@mm)
    .summarize(max_bytes: 8388608, batch_rows: 64)
print("Trace-flag days, then days with positive recorded precipitation")
print(trace.rows())
print(recorded_positive.rows())

// This is an illustrative threshold screen, not a heatwave definition or a
// climatological percentile baseline. Observation IDs remain exact YYYYMMDD
// values for this single fixed station; all native views retain source order.
let warm = accepted.filter(data.col("TMAX_DegC") >= 30.0@degC)
    .select(List("observation_id", "date", "TMAX_DegC", "TMIN_DegC", "PRCP_mm", "PRCP_mflag"))
    .collect(max_rows: 0, max_bytes: 8388608, batch_rows: 64)
print("Illustrative days with daily maximum at least 30 degC")
print(warm.rows())
print(warm.rows() == 39 && observations.rows() == 366)
print("Observed mean daily maximum and population spread")
print(basic.unit_mean("TMAX_DegC", "degC"))
print(basic.unit_std("TMAX_DegC", "delta_degC"))

// Daily dates are calendar labels, not invented midnight UTC instants. Day of
// year is an explicit plotting coordinate for this one leap year, not a new
// timestamp or a fitted time trend. All 366 rows are inspected in these charts.
let seasonal = plot.scatter(accepted, x: "day_of_year", y: "TMAX_DegC", max_points: 366)
    .title("NOAA Central Park, 2024: daily maximum temperature")
    .labels("Day of year (2024, leap year)", "Daily maximum temperature (degC)")
    .theme(plot.theme("light"))
    .size(1000, 600)
let shape = plot.distribution(deep, "TMAX_DegC", bins: 18)
    .title("NOAA Central Park, 2024: observed daily maxima")
    .theme(plot.theme("light"))
    .size(1000, 600)
let paired = plot.scatter(accepted, x: "TMIN_DegC", y: "TMAX_DegC", max_points: 366)
    .title("NOAA Central Park, 2024: paired daily temperatures")
    .theme(plot.theme("light"))
    .size(1000, 600)
let coverage = plot.missingness(deep)
    .title("NOAA Central Park, 2024: no missing selected values")
    .theme(plot.theme("light"))
    .size(1000, 600)
seasonal.show()
shape.show()
paired.show()
coverage.show()

// Persist these checkpoint values with the source CSV, source-months.dly,
// manifest, analysis program and build identity. stdout can be archived by the
// CLI; each chart also preserves structured lineage in its native export.
print("REPRODUCIBILITY CHECKPOINT: source, basic and deep lineage")
print(observations.provenance())
print(deep.profile_json())
print(warm.provenance())
print(seasonal.provenance())
print("Complete snapshot traversal does not erase truncated origin")
let preview = source.preview(7, 8388608, 64)
let inspected_prefix = preview.scan_table().profile(8388608, 64)
print(inspected_prefix.complete() && !inspected_prefix.origin_complete())
// One station-year is seasonal descriptive evidence. It does not establish a
// long-term climate trend, causal effect, heatwave risk or representative city
// exposure. Station history, observation windows, homogeneous long records and
// appropriate uncertainty models are necessary for those questions.