From 80c3d9d4f16076f04ea2da74724529126005f357 Mon Sep 17 00:00:00 2001 From: brodrick Date: Thu, 26 Feb 2026 11:43:30 -0800 Subject: [PATCH 1/2] add pca plots --- spectral_util/__main__.py | 1 + spectral_util/common/__init__.py | 10 ++++- spectral_util/common/plotting.py | 71 +++++++++++++++++++++++++++++++- 3 files changed, 80 insertions(+), 2 deletions(-) diff --git a/spectral_util/__main__.py b/spectral_util/__main__.py index 4325861..eee0099 100644 --- a/spectral_util/__main__.py +++ b/spectral_util/__main__.py @@ -15,6 +15,7 @@ def cli(): cli.add_command(mosaic.cli, name='mosaic') cli.add_command(common.cli_quicklook, name='quicklooks') cli.add_command(common.cli_plot, name='plot') +cli.add_command(common.cli_extract, name='extract') if __name__ == '__main__': cli() \ No newline at end of file diff --git a/spectral_util/common/__init__.py b/spectral_util/common/__init__.py index 1f1f728..48a3735 100644 --- a/spectral_util/common/__init__.py +++ b/spectral_util/common/__init__.py @@ -2,7 +2,8 @@ import click from .quicklooks import ndvi, rgb, nbr -from .plotting import plot_basic_overview +from .plotting import plot_basic_overview, plot_pcs +from .extract import extract_polygon_data, extract_polygon_stats, extract_polygon_to_geojson @click.group() def cli_quicklook(): @@ -12,7 +13,14 @@ def cli_quicklook(): def cli_plot(): pass +@click.group() +def cli_extract(): + pass + cli_quicklook.add_command(rgb) cli_quicklook.add_command(nbr) cli_quicklook.add_command(ndvi) cli_plot.add_command(plot_basic_overview) +cli_plot.add_command(plot_pcs) +# Extract commands are functions, not Click commands, so we don't add them here +# They are available as module-level functions via spectral_util.common.extract diff --git a/spectral_util/common/plotting.py b/spectral_util/common/plotting.py index 2a319be..24a1a45 100644 --- a/spectral_util/common/plotting.py +++ b/spectral_util/common/plotting.py @@ -7,6 +7,7 @@ from spectral_util.common.quicklooks import get_rgb import matplotlib.pyplot as plt from sklearn.cluster import MiniBatchKMeans +import scipy def plot_spectra(input_file, n_points=10, method='even', seed=13): """ @@ -71,7 +72,6 @@ def plot_spectra(input_file, n_points=10, method='even', seed=13): axes[1].set_xlabel("Wavelength") wl_nan = meta.wavelengths.copy() - print(spectra_to_plot[0,:]) wl_nan[np.isclose(spectra_to_plot[0,:], -0.01, atol=1e-5)] = np.nan for i, spec in enumerate(spectra_to_plot): axes[1].plot(wl_nan, spec, label=labels[i], c=cmap(i % cmap.N)) @@ -98,6 +98,75 @@ def plot_basic_overview(input_file, output_file, n_points, method): plt.savefig(output_file, bbox_inches='tight', dpi=300) click.echo(f"Plot saved to {output_file}") + +@click.command() +@click.argument('input_file', type=click.Path(exists=True)) +@click.option('--output_file', '-o', default=None, help='Path to save the plot (if not provided, will display instead)') +def plot_pcs(input_file, output_file): + """ + Visualizes an image and consistent spectra. + """ + # Click passes None if bands is not provided, handled in function + pc_figure(input_file) + if output_file: + plt.savefig(output_file, bbox_inches='tight', dpi=300) + click.echo(f"Plot saved to {output_file}") + + +def build_pca(x): + mu = x.mean(axis=0) + C = np.cov(x-mu, rowvar=False) + [v, d] = scipy.linalg.eigh(C) + + return d, mu + + +def do_pca(x, npca, d, mu, pca_offset=0): + # Project, redimension as an image with "npca" channels, and segment + if pca_offset == 0: + pca = (x - mu) @ d[:, -npca:] + else: + pca = (x - mu) @ d[:, -(npca+pca_offset):-1*pca_offset] + return pca[...,::-1] + + + +def pc_figure(input_file, start=0, stop=99, seed=13, n_points=10_000): + """ + Plots single pand pcs + + Args: + input_file (str): Path to the input spectral file. + start (int): Starting index of the principal components to plot. + stop (int): Ending index of the principal components to plot. + """ + np.random.seed(seed) + # Load data + meta, data = load_data(input_file) + build_data = data.copy().reshape((data.shape[0] * data.shape[1], data.shape[2])) + build_data = build_data[np.random.permutation(build_data.shape[0])[:n_points],:] + d, mu = build_pca(build_data) + pca_out = np.zeros((data.shape[0], data.shape[1], stop-start+1)) + n_pcs = stop - start + 1 + for _line in range(data.shape[0]): + pca_out[_line,...] = do_pca(data[_line,...], n_pcs, d, mu, pca_offset=start) + + plt.figure(figsize=(10, 10)) + rows = int(np.ceil(np.sqrt(n_pcs))) + cols = int(np.ceil(n_pcs / rows)) + grid = plt.GridSpec(rows, cols, wspace=0.13, hspace=0.13) + for i in range(n_pcs): + plt.subplot(grid[i]) + plt.imshow(pca_out[...,i], vmin=np.percentile(pca_out[...,i], 2), vmax=np.percentile(pca_out[...,i], 98)) + plt.title(f"PC {start+i}", fontsize=8) + plt.axis('off') + plt.tight_layout() + plt.show() + + + + + if __name__ == '__main__': plot_basic_overview() From 386feb48b0dc8bd042fb89ffe551044d7b33371b Mon Sep 17 00:00:00 2001 From: brodrick Date: Thu, 26 Feb 2026 11:51:45 -0800 Subject: [PATCH 2/2] remove lingering extract from separate pr --- spectral_util/common/__init__.py | 1 - 1 file changed, 1 deletion(-) diff --git a/spectral_util/common/__init__.py b/spectral_util/common/__init__.py index 48a3735..7ba334c 100644 --- a/spectral_util/common/__init__.py +++ b/spectral_util/common/__init__.py @@ -3,7 +3,6 @@ import click from .quicklooks import ndvi, rgb, nbr from .plotting import plot_basic_overview, plot_pcs -from .extract import extract_polygon_data, extract_polygon_stats, extract_polygon_to_geojson @click.group() def cli_quicklook():