import matplotlib matplotlib.use("Agg") # no display needed for plot smoke tests import numpy as np import pyarrow as pa import pyarrow.parquet as pq import pytest from giant.analysis import ( RAW_TARGET_NAMES, SampleCollection, constraint_report, constraint_report_pl, correlation_matrices, direction_alignment, load_predicted_local, marginal_table, marginal_table_pl, plot_constraint_violations, plot_correlation_matrices, plot_direction_alignment, plot_marginals, plot_pairwise, ) from giant.constants import ( LOCAL_TARGET_NAMES, PREDICT_COORD_METADATA_KEY, PREDICT_SCHEMA_VERSION, PREDICT_SCHEMA_VERSION_KEY, ) from giant.data.transforms import log_transform def _unit_vectors(rng, n): v = rng.standard_normal((n, 3)).astype(np.float32) return v / np.linalg.norm(v, axis=1, keepdims=True) def _make_collection(n=200, seed=0, gen_offset=0.0) -> SampleCollection: rng = np.random.default_rng(seed) real = np.column_stack( [ rng.uniform(0.1, 5.0, n), # step_length rng.uniform(0.1, 5.0, n), # delta_e rng.uniform(0.1, 5.0, n), # edep _unit_vectors(rng, n), # post_dir _unit_vectors(rng, n), # travel_dir ] ).astype(np.float32) gen = real + gen_offset pre_E = rng.uniform(1.0, 100.0, n).astype(np.float32) cond_cont_raw = np.column_stack( [ rng.standard_normal((n, 3)), pre_E, rng.standard_normal((n, 3)), rng.integers(0, 5, n), rng.integers(0, 3, n), ] ).astype(np.float32) return SampleCollection( cond_cont_raw=cond_cont_raw, pdg=rng.choice([11, -11, 22], size=n), material=rng.choice(["W", "Pb"], size=n), real_raw=real, gen_raw=gen, real_norm=real, gen_norm=gen, ) def test_marginal_table_aggregate_has_all_dims(): table = marginal_table(_make_collection()) assert set(table["dim"]) == set(RAW_TARGET_NAMES) assert (table["group"] == "all").all() @pytest.mark.parametrize("group_by", ["pdg", "material", "energy"]) def test_marginal_table_grouped_covers_all_rows(group_by): collection = _make_collection() table = marginal_table(collection, group_by=group_by) assert table["n"].groupby(table["group"]).first().sum() == len(collection.pdg) def test_marginal_table_identical_distributions_have_zero_kl(): collection = _make_collection(gen_offset=0.0) table = marginal_table(collection) np.testing.assert_allclose(table["kl_real_gen"], 0.0, atol=1e-6) def test_marginal_table_shifted_distribution_has_positive_kl(): collection = _make_collection(gen_offset=3.0) table = marginal_table(collection) assert (table["kl_real_gen"] > 0).all() def test_correlation_matrices_are_symmetric_unit_diagonal(): real_corr, gen_corr = correlation_matrices(_make_collection()) for corr in (real_corr, gen_corr): np.testing.assert_allclose(np.diag(corr), 1.0, atol=1e-5) np.testing.assert_allclose(corr, corr.T, atol=1e-5) def test_direction_alignment_real_data_is_unit_norm_dot_product(): real_cos, gen_cos = direction_alignment(_make_collection()) assert np.all(real_cos >= -1.0 - 1e-5) and np.all(real_cos <= 1.0 + 1e-5) assert np.all(gen_cos >= -1.0 - 1e-5) and np.all(gen_cos <= 1.0 + 1e-5) def test_constraint_report_clean_data_has_no_violations(): report = constraint_report(_make_collection(gen_offset=0.0)) assert (report["violation_rate"] == 0.0).all() def test_constraint_report_flags_negative_log_dims_and_bad_norms(): collection = _make_collection(gen_offset=0.0) collection.gen_raw[:, 0] = -1.0 # negative step_length collection.gen_raw[:, 3:6] *= 2.0 # post_dir no longer unit norm report = constraint_report(collection) violations = dict(zip(report["check"], report["violation_rate"])) assert violations["step_length >= 0"] == 1.0 assert violations["post_dir unit norm"] == 1.0 def test_plot_marginals_runs_without_error(): fig = plot_marginals(_make_collection()) assert fig is not None def test_plot_marginals_grouped_runs_without_error(): fig = plot_marginals(_make_collection(), group_by="material") assert fig is not None def test_plot_correlation_matrices_runs_without_error(): fig = plot_correlation_matrices(_make_collection()) assert fig is not None def test_plot_pairwise_runs_without_error(): fig = plot_pairwise(_make_collection()) assert fig is not None def test_plot_direction_alignment_runs_without_error(): fig = plot_direction_alignment(_make_collection()) assert fig is not None def test_plot_constraint_violations_runs_without_error(): fig = plot_constraint_violations(_make_collection()) assert fig is not None def _write_predicted_local_parquet(path, n=50, metadata=None, rng=None): """Mimic `giant predict --coord local`'s output schema for the loader tests.""" rng = rng or np.random.default_rng(0) true_log_local = rng.standard_normal((n, 9)).astype(np.float32) true_log_local[:, :3] = log_transform( rng.uniform(0.1, 5.0, (n, 3)).astype(np.float32) ) pred_log_local = true_log_local + rng.normal(0, 0.01, (n, 9)).astype(np.float32) table = pa.table( { "event_id": rng.integers(0, 10, n), "pdg": rng.choice([11, -11, 22], n), "pre_x": rng.standard_normal(n).astype(np.float32), "pre_y": rng.standard_normal(n).astype(np.float32), "pre_z": rng.standard_normal(n).astype(np.float32), "pre_E": rng.uniform(1.0, 100.0, n).astype(np.float32), "pre_dx": rng.standard_normal(n).astype(np.float32), "pre_dy": rng.standard_normal(n).astype(np.float32), "pre_dz": rng.standard_normal(n).astype(np.float32), "material": rng.choice(["W", "Pb"], n), "layer_id": rng.integers(0, 10, n).astype(np.int32), "n_sec": rng.integers(0, 3, n).astype(np.int32), **{ f"pred_{name}": pred_log_local[:, j] for j, name in enumerate(LOCAL_TARGET_NAMES) }, **{ f"true_{name}": true_log_local[:, j] for j, name in enumerate(LOCAL_TARGET_NAMES) }, } ) if metadata is not None: table = table.replace_schema_metadata(metadata) pq.write_table(table, path) return true_log_local, pred_log_local def test_load_predicted_local_round_trips_values(tmp_path): path = tmp_path / "predicted_local.parquet" true_log_local, pred_log_local = _write_predicted_local_parquet( path, metadata={ PREDICT_COORD_METADATA_KEY: "local", PREDICT_SCHEMA_VERSION_KEY: PREDICT_SCHEMA_VERSION, }, ) collection = load_predicted_local(path) expected_real = true_log_local.copy() expected_real[:, :3] = np.exp(expected_real[:, :3]) - 1e-8 expected_gen = pred_log_local.copy() expected_gen[:, :3] = np.exp(expected_gen[:, :3]) - 1e-8 np.testing.assert_allclose(collection.real_raw, expected_real, atol=1e-4) np.testing.assert_allclose(collection.gen_raw, expected_gen, atol=1e-4) assert collection.real_norm is None assert collection.gen_norm is None def test_load_predicted_local_usable_by_downstream_plots(tmp_path): path = tmp_path / "predicted_local.parquet" _write_predicted_local_parquet( path, metadata={ PREDICT_COORD_METADATA_KEY: "local", PREDICT_SCHEMA_VERSION_KEY: PREDICT_SCHEMA_VERSION, }, ) collection = load_predicted_local(path) assert marginal_table(collection) is not None assert plot_marginals(collection) is not None def test_load_predicted_local_rejects_missing_metadata(tmp_path): path = tmp_path / "no_metadata.parquet" _write_predicted_local_parquet(path, metadata=None) with pytest.raises(ValueError, match="no '.*' parquet metadata"): load_predicted_local(path) def test_load_predicted_local_rejects_global_coord(tmp_path): path = tmp_path / "global.parquet" _write_predicted_local_parquet( path, metadata={ PREDICT_COORD_METADATA_KEY: "global", PREDICT_SCHEMA_VERSION_KEY: PREDICT_SCHEMA_VERSION, }, ) with pytest.raises(ValueError, match="coord=local"): load_predicted_local(path) def test_load_predicted_local_rejects_mismatched_schema_version(tmp_path): path = tmp_path / "old_version.parquet" _write_predicted_local_parquet( path, metadata={ PREDICT_COORD_METADATA_KEY: "local", PREDICT_SCHEMA_VERSION_KEY: "999", }, ) with pytest.raises(ValueError, match="schema version"): load_predicted_local(path) def _predicted_local_path(tmp_path, n=200, seed=0): path = tmp_path / "predicted_local.parquet" _write_predicted_local_parquet( path, n=n, rng=np.random.default_rng(seed), metadata={ PREDICT_COORD_METADATA_KEY: "local", PREDICT_SCHEMA_VERSION_KEY: PREDICT_SCHEMA_VERSION, }, ) return path @pytest.mark.parametrize("group_by", [None, "pdg", "material", "energy"]) def test_marginal_table_pl_matches_numpy_version(tmp_path, group_by): path = _predicted_local_path(tmp_path) collection = load_predicted_local(path) expected = marginal_table(collection, group_by=group_by).sort_values( ["group", "dim"] ) actual = ( marginal_table_pl(path, group_by=group_by).sort(["group", "dim"]).to_pandas() ) assert list(expected["group"]) == list(actual["group"]) assert list(expected["n"]) == list(actual["n"]) for col in ["real_mean", "gen_mean", "real_std", "gen_std", "kl_real_gen"]: np.testing.assert_allclose( expected[col].to_numpy(), actual[col].to_numpy(), atol=1e-4, rtol=1e-4, ) def test_marginal_table_pl_rejects_missing_metadata(tmp_path): path = tmp_path / "no_metadata.parquet" _write_predicted_local_parquet(path, metadata=None) with pytest.raises(ValueError, match="no '.*' parquet metadata"): marginal_table_pl(path) def test_constraint_report_pl_matches_numpy_version(tmp_path): path = _predicted_local_path(tmp_path) collection = load_predicted_local(path) expected = constraint_report(collection) actual = constraint_report_pl(path).to_pandas() assert list(expected["check"]) == list(actual["check"]) np.testing.assert_allclose( expected["violation_rate"].to_numpy(), actual["violation_rate"].to_numpy(), atol=1e-6, ) np.testing.assert_allclose( expected["mean_abs_error"].to_numpy(), actual["mean_abs_error"].to_numpy(), atol=1e-4, )