{ "cells": [ { "cell_type": "markdown", "id": "d0fe5552-0348-4e79-bbf7-fe37b7edadda", "metadata": {}, "source": [ "# Notebook replicating a typical hyperspectral datacube analysis done using pyfresco" ] }, { "cell_type": "code", "execution_count": 1, "id": "cff62ae2-346a-4cd9-a9b1-dd2b4ff01526", "metadata": {}, "outputs": [], "source": [ "import pyfresco" ] }, { "cell_type": "markdown", "id": "b9d345e4-99bb-4f29-8efc-c52a790d2f44", "metadata": {}, "source": [ "## Data import\n", "\n", "In the following code cell we define the paths pointing to the img and header files of the proper hyperspectral datacube and of the summary parameters datacube. These two datacubes are the ones needed by pyFRESCO to work properly." ] }, { "cell_type": "code", "execution_count": 2, "id": "210b3f99-3280-47ef-abd1-c4c0756b3381", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "['R770', 'RBR', 'BD530_2', 'SH600_2', 'SH770', 'BD640_2', 'BD860_2', 'BD920_2', 'RPEAK1', 'BDI1000VIS', 'R440', 'IRR1', 'BDI1000IR', 'OLINDEX3', 'R1330', 'BD1300', 'LCPINDEX2', 'HCPINDEX2', 'VAR', 'ISLOPE1', 'BD1400', 'BD1435', 'BD1500_2', 'ICER1_2', 'BD1750_2', 'BD1900_2', 'BD1900R2', 'BDI2000', 'BD2100_2', 'BD2165', 'BD2190', 'MIN2200', 'BD2210_2', 'D2200', 'BD2230', 'BD2250', 'MIN2250', 'BD2265', 'BD2290', 'D2300', 'BD2355', 'SINDEX2', 'ICER2_2', 'MIN2295_2480', 'MIN2345_2537', 'BD2500_2', 'BD3000', 'BD3100', 'BD3200', 'BD3400_2', 'CINDEX2', 'BD2600', 'IRR2', 'IRR3', 'R530', 'R600', 'R1080', 'R1506', 'R2529', 'R3920']\n" ] } ], "source": [ "path_if_mtrdr = \"FRT00009b5a/frt00009b5a_07_if165j_mtr3.img\" # Path to the proper hyperspectral datacube\n", "head_if_mtrdr = \"FRT00009b5a/frt00009b5a_07_if165j_mtr3.hdr\" # Path to the header of the proper hyperspectral datacube\n", "\n", "path_sr_mtrdr = \"FRT00009b5a/frt00009b5a_07_sr165j_mtr3.img\" # Path to the summary parameters datacube\n", "head_sr_mtrdr = \"FRT00009b5a/frt00009b5a_07_sr165j_mtr3.hdr\" # Path to the header of the summary parameters datacube\n", "\n", "img , img_sr , wavelength , sr_names = pyfresco.open_raw(path_if_mtrdr , head_if_mtrdr , path_sr_mtrdr , head_sr_mtrdr) # Open the two datacubes\n", "\n", "Nbands = len(wavelength)\n", "print(sr_names)" ] }, { "cell_type": "markdown", "id": "db2e30f5-dcaa-49ee-be5e-2a4d68630af8", "metadata": {}, "source": [ "## Creation of the false color RGB map\n", "\n", "In the two following code cells we use pyFRESCO to create the false color and phyllosilicates RGB maps with methodology defined in Viviano-Beck et al., 2014." ] }, { "cell_type": "code", "execution_count": 3, "id": "2233e240-bc10-4b33-a2d2-65b75d3bc55d", "metadata": { "scrolled": true }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[58, 57, 56] ['R2529', 'R1506', 'R1080']\n", "[58, 57, 56] ['R2529', 'R1506', 'R1080']\n", "[58, 57, 56] ['R2529', 'R1506', 'R1080']\n", "Red channel limits = [ 0.145 , 0.3 ]\n", "Green channel limits = [ 0.16 , 0.3 ]\n", "Blue channel limits = [ 0.15 , 0.28 ]\n" ] } ], "source": [ "FAL = pyfresco.RGBImageManipulator(img, img_sr, 'FAL', 0, 0, 0, wavelength) # Define the object FAL\n", "\n", "FALmap, selected_stretch = FAL.RGBmapmake(FALSE = FAL, bi = 100, clip = True, cumhist = False,\n", " preset_true_colors = False, use_false_color = False) # Effective creation of the false color RGB map\n", "\n", "FAL.savemap('FAL' , folder = 'FRT00009b5a' , extension = '.tif' , show = True) # saving of the RGB map in a given folder, both in txt and tif files\n", "\n", "print('Red channel limits = ', '[', selected_stretch[0] , ',' , selected_stretch[3] , ']')\n", "print('Green channel limits = ', '[', selected_stretch[1] , ',' , selected_stretch[4] , ']')\n", "print('Blue channel limits = ', '[', selected_stretch[2] , ',' , selected_stretch[5] , ']')" ] }, { "cell_type": "markdown", "id": "a029f2c2-45d2-47f1-a8e3-b2a778e421db", "metadata": {}, "source": [ "The output of the above cell should look like this: ![False color map](notebook1_images/fal_making.png)" ] }, { "cell_type": "code", "execution_count": 6, "id": "6943020c-5d73-4398-b7d2-40ca5ed37bec", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[39, 33, 26] ['D2300', 'D2200', 'BD1900R2']\n", "[39, 33, 26] ['D2300', 'D2200', 'BD1900R2']\n", "Red channel limits = [ 0 , 0.0105 ]\n", "Green channel limits = [ 0 , 0.008 ]\n", "Blue channel limits = [ 0.005 , 0.027 ]\n" ] } ], "source": [ "PHY = pyfresco.RGBImageManipulator(img, img_sr, 'PHY', 0, 0, 0, wavelength) # Definition of the object PHY\n", "\n", "phy, selected_stretch = PHY.RGBmapmake(FALSE = FALmap, bi = 500, # Effective creation of the phyllosilicate RGB map\n", " clip = True, cumhist = False,\n", " preset_true_colors = False, \n", " use_false_color = True,\n", " R_min_in=[0,0.3] , R_max_in=[0,0.3] , \n", " G_min_in=[0,0.03] , G_max_in=[0,0.03] , \n", " B_min_in=[0,0.03] , B_max_in=[0,0.03] , \n", " slider_step=0.0005)\n", "\n", "PHY.savemap('PHY' , folder = 'FRT00009b5a' , extension = '.tif' , show = True) # Saving of the RGB map ina given folder as in the cell above\n", "\n", "print('Red channel limits = ', '[', selected_stretch[0] , ',' , selected_stretch[3] , ']')\n", "print('Green channel limits = ', '[', selected_stretch[1] , ',' , selected_stretch[4] , ']')\n", "print('Blue channel limits = ', '[', selected_stretch[2] , ',' , selected_stretch[5] , ']')" ] }, { "cell_type": "markdown", "id": "1da1a63e-93aa-44a2-97b2-30583669aaee", "metadata": {}, "source": [ "The output of the above cell should look like this: ![Phyllosilicates RGB map](notebook1_images/phy_making.png)" ] }, { "cell_type": "markdown", "id": "88a9e818-95a4-441d-b19d-de4ff0803f35", "metadata": {}, "source": [ "## Extraction of spectra from the datacube\n", "\n", "In the following cells we extract a spectrum from the hyperspectral datacube using the phyllosilicate RGB map as a guideline.\n", "\n", "The spectrum will be the mean spectrum of a set of spectra inside a Region of Interest (ROI), manually selected by drawing a polygon on the phyllosilicate RGB map." ] }, { "cell_type": "code", "execution_count": 7, "id": "b18130a1-d923-4608-8c1f-4f8b72e6f1a8", "metadata": {}, "outputs": [], "source": [ "SE = pyfresco.SpectraExtract(img , Nbands , wavelength , 750 , 2600) # Definition of the object for spectral extraction\n", "hyd = SE.upload_map('PHY' , folder = 'FRT00009b5a') # Upload an RGB map for the ROI selection" ] }, { "cell_type": "code", "execution_count": 8, "id": "d65bd7af-adf9-4123-8e27-9919c4ef816e", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Number of points inside the drewn polygon: 672\n" ] } ], "source": [ "polygon_spectra , L , mask = SE.polygon_spectra(save_pixel = True , # Select polygonal ROI on the RGB map\n", " folder = 'FRT00009b5a' ,\n", " name = 'phy_ROI_coordinates')\n", "polygon_spectrum , polygon_error = SE.final_spectra() # Compute mean and standard deviation spectrum of the ROI and plot it" ] }, { "cell_type": "markdown", "id": "0db1da15-0e7d-420a-8249-7ab0ac4d72d9", "metadata": {}, "source": [ "The outputs of the above cell should something like the following images:\n", "\n", "Selection of the ROI by polygon:\n", "![ROI selection](notebook1_images/phy_ROI.png)\n", "\n", "Resulting mean spectrum with stadard deviation:\n", "![Mean spectrum of the ROI](notebook1_images/ROI_mean.png)" ] }, { "cell_type": "markdown", "id": "7e0748cc-5b39-44c5-93fa-b6ee787b0770", "metadata": {}, "source": [ "## Normalization of target ROI spectrum\n", "\n", "PyFRESCO offers various ways for the normalization of the spectrum. Following we use a 'classically' used method that we can call 'Neutral Polygon Method', in which a median spectrum from a dark ROI is selected to then perform the following operation to obtain a normalized (or ratioed) spectrum:\n", "\n", "$$ S_{Normalized} = \\frac{S_{ROI}}{S_{Neutral}} $$\n", "\n", "The further propagated error can be obtained either by classic error propagation or via a bootstrap method." ] }, { "cell_type": "code", "execution_count": 9, "id": "73748ee4-9f2e-4f7d-943f-5560438c8f07", "metadata": {}, "outputs": [], "source": [ "SN = pyfresco.SpectraNorm(hyd , Nbands , img , img_sr , wavelength , # Definition of the object for spectral normalization\n", " polygon_spectrum , polygon_error , polygon_spectra ,\n", " 750 , 2600)" ] }, { "cell_type": "code", "execution_count": 10, "id": "0ef3393a-067f-4e0e-b344-a72eb3d54516", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Number of points inside the drewn polygon: 1450\n", "(672, 489) (1450, 489)\n", "0.7220216606498195 % of the errors, set to zero for simplicity, are in reality NaN values.\n" ] } ], "source": [ "n_polygon_spectra , n_polygon , n_polygon_err , Ln , mask = SN.neutral_polygon_spectra() # Extraction of a polygon\n", "\n", "print(polygon_spectra.shape , n_polygon_spectra.shape)\n", "\n", "SN.plot_together() # Plot of target and neutral spectra\n", "normpoly , errorpoly = SN.bootstrapnorm(convexhull = False , N = 5000) # Computation of the normalized spectrum with bootstrap error\n", "SN.normplot() # Plot fo the normalized spectrum\n", " \n", "polynorm , polyerr = SN.norm_spectra() # Computation of the classically normalized spectrum\n", "SN.normplot() # Plot of the normalized spectrum" ] }, { "cell_type": "markdown", "id": "cdae81d0-b20d-4638-9134-4768943eb6c6", "metadata": {}, "source": [ "The outputs of the above cell should be something like:\n", "\n", "Selection of the neutral ROI with a polygon:\n", "![Selection of the neutral ROI](notebook1_images/neutral_ROI.png)\n", "\n", "Target ROI mean spectrum and neutral ROI median spectrum:\n", "![Target and neutral spectra together](notebook1_images/neutral_and_mean.png)\n", "\n", "Resulting normalized spectrum with bootstrap derived error:\n", "![Bootstrap normalized spectrum](notebook1_images/bootstrap_norm.png)\n", "\n", "Resulting normalized spectrum with classically derived error:\n", "![Classically normalized spectrum](notebook1_images/classic_norm.png)" ] }, { "cell_type": "markdown", "id": "d0ba1434-8bea-4a12-89ee-5e8091543a9a", "metadata": {}, "source": [ "## Analysis of the normalized spectrum\n", "\n", "PyFRESCO offers various ways to ntaively analyze the normalized spectrum. Following, we use a comparison method that is based on the comparison of the absorption band positions with other, validated spectra taken from the MICA files (Viviano-Beck et al., 2015). \n", "\n", "The comparison can be done in three ways, from less complete to more complete ones, of which the last one offers also the possibility of automatically inferring line positions given a MICA files spectrum to compare." ] }, { "cell_type": "code", "execution_count": 11, "id": "bdd3fe23-2819-4235-a2a1-a8f05ee9afad", "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "C:\\Users\\mbaro\\miniconda3\\envs\\pyfresco\\lib\\site-packages\\pyfresco\\SpectraAnalysis.py:257: ParserWarning: Falling back to the 'python' engine because the 'c' engine does not support sep=None with delim_whitespace=False; you can avoid this warning by specifying engine='python'.\n", " spectra_info = pd.read_csv('spectra_MICA_LAB_info.csv' , sep = None ,\n" ] }, { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
Mineral NameMineral Typetxt name CRISMtxt name LABCRISM spectra cubeLab. MineralType of Lab. From MICAMICA Lab.CodeType of Lab. From websiteNotable Absorptions CRISM [nm]Notable Absorptions LAB [nm]Sample Grain Size [mum]
13Mg SmectitePhyllosilicatecrism_spec_mg_smectite.txtmg_smectite_LAB.txtFRT00009365SaponiteRELABLASA58GEO[1410 , 1920 , 2310 , 2390][1390 , 1910 , 2320 , 2390]max 4
\n", "
" ], "text/plain": [ " Mineral Name Mineral Type txt name CRISM \\\n", "13 Mg Smectite Phyllosilicate crism_spec_mg_smectite.txt \n", "\n", " txt name LAB CRISM spectra cube Lab. Mineral \\\n", "13 mg_smectite_LAB.txt FRT00009365 Saponite \n", "\n", " Type of Lab. From MICA MICA Lab.Code Type of Lab. From website \\\n", "13 RELAB LASA58 GEO \n", "\n", " Notable Absorptions CRISM [nm] Notable Absorptions LAB [nm] \\\n", "13 [1410 , 1920 , 2310 , 2390] [1390 , 1910 , 2320 , 2390] \n", "\n", " Sample Grain Size [mum] \n", "13 max 4 " ] }, "execution_count": 11, "metadata": {}, "output_type": "execute_result" } ], "source": [ "SA = pyfresco.SpectraAnalysis(normpoly , errorpoly , # Definition of the object for spectral analysis\n", " m_spec = polygon_spectrum , err_spec = polygon_error , \n", " n_spectra = n_polygon , n_err = n_polygon_err , wavelength = wavelength , \n", " MIN = 750 , MAX = 2600 , folder = None)\n", "\n", "SA.dataset_interaction('Mg Smectite' , use = 'Both') # Displaying details of our compositional guess" ] }, { "cell_type": "code", "execution_count": 12, "id": "52bd5316-7afd-4113-b4c3-30a7b6a3e255", "metadata": {}, "outputs": [], "source": [ "SA.simple_compare('Mg Smectite', use = 'Both') # Plot of the normalized spectrum with superimposed tabulated and validated absorption lines" ] }, { "cell_type": "markdown", "id": "4e17d9bd-2e7d-462a-a7bf-005092c1b7ce", "metadata": {}, "source": [ "The result of the above cell will be something like:\n", "\n", "![Simple comparison](notebook1_images/simple_comparison.png)" ] }, { "cell_type": "code", "execution_count": 13, "id": "03dd9135-3d18-4871-91c2-bcd15550cc3f", "metadata": {}, "outputs": [], "source": [ "SA.compare_spectra('Mg Smectite' , use = 'Both') # Plot of the normalized spectrum with superimposed tabulated and validated absorption lines along \n", " # with the CRISM and laboratory spectra from the MICA files" ] }, { "cell_type": "markdown", "id": "003af3e5-414b-4a0c-90e7-7637c93ea16f", "metadata": {}, "source": [ "The result of the above cell will be something like:\n", "\n", "![Spectra comparison](notebook1_images/compare_spectra.png)" ] }, { "cell_type": "code", "execution_count": 14, "id": "8f35e1f9-0d96-43f5-8dee-e47daa315f1b", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "CRISM derived absorption(s) = [np.float64(1414.59), np.float64(1914.87), np.float64(2311.18), np.float64(2390.58)]\n", "CRISM derived absorption(s) difference = [np.float64(4.589999999999918), np.float64(5.130000000000109), np.float64(1.1799999999998363), np.float64(0.5799999999999272)]\n", "LAB derived absorption(s) =. [np.float64(1414.59), np.float64(1914.87), np.float64(2311.18), np.float64(2390.58)]\n", "LAB derived absorption(s) difference = [np.float64(24.589999999999918), np.float64(4.869999999999891), np.float64(8.820000000000164), np.float64(0.5799999999999272)]\n" ] } ], "source": [ "absoC , absoL , diffC , diffL , CRISM_abs , LAB_abs = SA.compare_lines('Mg Smectite' , ran = 5 , s = 1 , t = 1) # Automatical inferring\n", " # of absorption line positions\n", "\n", "print('CRISM derived absorption(s) = ' , absoC)\n", "print('CRISM derived absorption(s) difference = ' , diffC)\n", "print('LAB derived absorption(s) =. ' , absoL)\n", "print('LAB derived absorption(s) difference = ' , diffL)" ] }, { "cell_type": "code", "execution_count": 15, "id": "885c5fca-1c03-46e8-baa3-d8a782dbc1ef", "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "C:\\Users\\mbaro\\miniconda3\\envs\\pyfresco\\lib\\site-packages\\pyfresco\\SpectraAnalysis.py:699: ParserWarning: Falling back to the 'python' engine because the 'c' engine does not support sep=None with delim_whitespace=False; you can avoid this warning by specifying engine='python'.\n", " spectra_info = pd.read_csv('spectra_MICA_LAB_info.csv' , sep = None ,\n" ] } ], "source": [ "SA.compare2(0 , n_polygon , n_polygon_err , normpoly , 'Mg Smectite') # Plot of the comaprison using the inferred lines" ] }, { "cell_type": "markdown", "id": "ccc1190d-838f-40f0-853a-dcffd3dfc2c3", "metadata": {}, "source": [ "The result of the above cells will be something like:\n", "\n", "![Complete comparison](notebook1_images/complete_comparison.png)" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (ipykernel)", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.10.20" } }, "nbformat": 4, "nbformat_minor": 5 }