From f9f3ad0712593b46bb81ba3c6dedadfd82b48e0d Mon Sep 17 00:00:00 2001 From: Abdussalam SALEH MGHIR <salehma@visu01.sis.cnes.fr> Date: Mon, 27 Sep 2021 07:42:43 +0000 Subject: [PATCH] ENH : Python script to compute average spectra in range and in azimuth for a Tiff image --- python_src/utils/compute_spectra_ImageTiff.py | 109 ++++++++++++++++++ 1 file changed, 109 insertions(+) create mode 100755 python_src/utils/compute_spectra_ImageTiff.py diff --git a/python_src/utils/compute_spectra_ImageTiff.py b/python_src/utils/compute_spectra_ImageTiff.py new file mode 100755 index 0000000..fb04409 --- /dev/null +++ b/python_src/utils/compute_spectra_ImageTiff.py @@ -0,0 +1,109 @@ +# coding:utf-8 +import sys +import numpy as np +import matplotlib.pyplot as plt +import os +import argparse +try: + import gdal +except ImportError: + import osgeo.gdal as gdal +def calcul_spectres_tiff(ficima,output_path, output_filename_az, output_filename_dist): + """ + Compute the average spectra in range and in azimuth from a complex SAR image (mandatory format .tif) + """ + # Average spectrum in azimut + ################################### + + print("Compute the average spectrum in azimut after loading the modules python and gdal") + #Read the image tiff + zone = read_tiff(ficima) + print ("Dimensions of the tif file ",np.shape(zone)) + (nb_lig, nb_col) = np.shape(zone[:,:,0]) + print ("nb lig = ",nb_lig," nb col = ", nb_col) + + lig = np.zeros(nb_lig, dtype = complex) + spectre_az = np.zeros(nb_lig, dtype = 'float') + + for i in range(nb_col): + lig = zone[:,i,0] + 1j*zone[:,i,1] #All the lines in complex + spectre_az = spectre_az + np.abs(np.fft.fft(lig))**2 + + spectre_az = np.sqrt(spectre_az) / float(nb_col) + + #Display and Save the spectrum + display_spectrum(spectre_az, "Average spectrum in azimut", output_path, output_filename_az) + + # Compute the average spectrum in range + ##################################### + + print("Compute the average spectrum in range") + col = np.zeros(nb_col, dtype = complex) + spectre_dist = np.zeros(nb_col, dtype = 'float') + + for i in range(nb_lig): + col = zone[i,:,0] + 1j*zone[i,:,1] #All the columns in complex + spectre_dist = spectre_dist + np.abs(np.fft.fft(col))**2 + + spectre_dist = np.sqrt(spectre_dist) / float(nb_lig) + + #Display and Save the spectrum + print("Display and Save the spectrum") + display_spectrum(spectre_dist, "Average spectrum in range", output_path, output_filename_dist) + + plt.show() + +def read_tiff(filename): + """ + Create a numpy array from the input image tiff. + """ + ds = gdal.Open(filename) + + nb_lig = ds.RasterYSize + nb_col = ds.RasterXSize + nb_bands = ds.RasterCount + + zone = np.zeros((nb_lig, nb_col, nb_bands)) + + for i in range(nb_bands): + zone[:,:,i] = np.array(ds.GetRasterBand(i+1).ReadAsArray()) + + return zone + +def display_spectrum(spectre, title, path, filename): + """ + Display the normalized spectrum. + """ + nb_val = np.size(spectre) + xvals = np.arange(nb_val) / float(nb_val) + + plt.figure() + plt.plot(xvals, spectre) + plt.ylim(0, 1.1*max(spectre)) + plt.title(title + "\n" + filename) + plt.savefig(os.path.join(path,filename)) + plt.draw() + + +if __name__ == "__main__": + description = "Compute of the spectra of the image tif" + parser = argparse.ArgumentParser(description=description) + parser.add_argument("path_input", type=str, help="Path of the input Tiff Image to get its spectrum. There is not its extension in the path") + parser.add_argument("path_output_dir", type=str, help="Path of the output directory where the spectra of the input will be saved") + args = parser.parse_args() + + + print("Usage : python compute_spectra_ImageTiff.py <path_input_tif_image without the extension file> <path_output_directory>") + print("The Tiff Image must have 2 bands in float/int. The script cannot read the image with one band in complex") + path_input = args.path_input + file_input = path_input.split('/')[-1] + input_path = os.path.dirname(path_input) + if(input_path==''): + input_path = '.' + nomima = file_input.split('.')[0] + ficima=os.path.join(input_path,nomima+".tif") + output_path=args.path_output_dir + output_filename_az="spectrum_az_"+nomima + output_filename_dist="spectrum_dist_"+nomima + + calcul_spectres_tiff(ficima,output_path, output_filename_az, output_filename_dist) -- GitLab