# -*- coding: utf-8 -*- # class for Jamet Neural Network Kd inversion. # Python version : # 2016-05-11 D.Dessailly david.dessailly@univ-littoral.fr # LOG (Laboratoire d'Oceanoligie et Geoscience) import numpy as np import sys import os Sensor = {'idSEAWIFS' : 0, 'idMODIS' : 1, 'idMERIS' : 2} NBWL = 6 # # Class to handle Jamet's Kd processing for MERIS/OLCI, SEAWIFS and MODIS Sensor # class KdJamet: # -- default __init__ Method +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ # Input : # idSensor : valid values -> ['idSEAWIFS', 'idMODIS', 'idMERIS'] # -- def __init__(self, idSensor): self.sensor = idSensor try: self.lutpath = os.getenv('KD_LUTS_PATH') except: print "KD_LUTS_PATH not set, add it to your environment variables\n for bash user add the following line to your .bashrc file :\n export KD_LUTS_PATH=/where/are/my/KD/LUTS\n then:\n source ~/.bashrc" sys.exit() if idSensor in Sensor: # self.init_variables() try: self.read_LUTs_Kd() except: print "Failed to read %s LUTs" % (idSensor[2:]) raise sys.exit() else: print "invalid Sensor ID : %s\n valid values : idSEAWIFS, idMODIS, idMERIS" % (idSensor) sys.exit() # -- METHOD init_variables ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ # Sensor dependants variables and Perceptron parameters setting # -- def init_variables(self): print "class KdJamet.init_variables Method" # ======================================================================================================= # MERIS / OLCI if self.sensor == 'idMERIS': self.WL = np.array([413, 443, 490, 510, 560, 665]) # Variables definiton for network Rrs490/Rrs555 >= .85 ------------------------------------------------- # hidden layer number self.rsup_NC=2 # input data number self.rsup_NE=5 # number of neuron in the first hidden layer self.rsup_NC1=7 # number of neuron in the second hidden layer self.rsup_NC2=4 # nb output self.rsup_NS=1 # nb input + nb output*/ self.rsup_NES=6 # LUTs file names self.rsup_LUT_POIDS = "KdOLCI_poids_IOCCG_NOMAD_BOUM_ss_12_COASTCOLOUR_443_to_555_log_Kd_lambda_merge_seuil_15_angle_ascii_6x1_hh_7_4_switch_110416.sn" self.rsup_LUT_MOY = "Moy_KdOLCI_IOCCG_NOMAD_BOUM_ss_12_COASTCOLOUR_412_to_670_log_Kd_lambda_merge_seuil_15_angle_publi_ss_620.dat" self.rsup_LUT_ECART = "Ecart_KdOLCI_IOCCG_NOMAD_BOUM_ss_12_COASTCOLOUR_412_to_670_log_Kd_lambda_merge_seuil_15_angle_publi_ss_620.dat" # Variables definiton for network < .85 ----------------------------------- # nombre de couches cachees - hidden layer number self.rinf_NC=2 # input data number*/ self.rinf_NE=6 # number of neuron in the first hidden layer self.rinf_NC1=5 # number of neuron in the second hidden layer self.rinf_NC2=5 # nb output self.rinf_NS=1 # nb input + nb output self.rinf_NES=7 #LUTs file names self.rinf_LUT_POIDS = "KdOLCI_poids_IOCCG_NOMAD_BOUM_ss_12_COASTCOLOUR_443_to_670_log_Kd_lambda_merge_seuil_15_angle_ascii_6x1_hh_5_5_switch_110416.sn" self.rinf_LUT_MOY = "Moy_KdOLCI_IOCCG_NOMAD_BOUM_ss_12_COASTCOLOUR_412_to_670_log_Kd_lambda_merge_seuil_15_angle_publi_ss_620.dat" self.rinf_LUT_ECART = "Ecart_KdOLCI_IOCCG_NOMAD_BOUM_ss_12_COASTCOLOUR_412_to_670_log_Kd_lambda_merge_seuil_15_angle_publi_ss_620.dat" # ======================================================================================================= # MODIS if self.sensor == 'idMODIS': self.WL = np.array([412, 443, 488, 531, 547, 667]) # Variables definiton for network Rrs490/Rrs555 >= .85 ------------------------------------------------- # hidden layer number self.rsup_NC=2 # input data number self.rsup_NE=5 # number of neuron in the first hidden layer self.rsup_NC1=5 # number of neuron in the second hidden layer self.rsup_NC2=4 # nb output self.rsup_NS=1 # nb input + nb output*/ self.rsup_NES=6 # LUTs file names self.rsup_LUT_POIDS = "KdMODIS_poids_IOCCG_NOMAD_BOUM_ss_12_COASTCOLOUR_443_to_555_log_Kd_lambda_merge_seuil_15_angle_ascii_6x1_hh_5_4_publi_CSIRO.sn" self.rsup_LUT_MOY = "Moy_KdMODIS_IOCCG_NOMAD_BOUM_ss_12_COASTCOLOUR_412_to_555_log_Kd_lambda_merge_seuil_15_ss_645_publi_sup_CSIRO.dat" self.rsup_LUT_ECART = "Ecart_KdMODIS_IOCCG_NOMAD_BOUM_ss_12_COASTCOLOUR_412_to_555_log_Kd_lambda_merge_seuil_15_ss_645_publi_sup_CSIRO.dat" # Variables definiton for network < .85 ----------------------------------- # nombre de couches cachees - hidden layer number self.rinf_NC=2 # input data number*/ self.rinf_NE=6 # number of neuron in the first hidden layer self.rinf_NC1=10 # number of neuron in the second hidden layer self.rinf_NC2=7 # nb output self.rinf_NS=1 # nb input + nb output self.rinf_NES=7 #LUTs file names self.rinf_LUT_POIDS = "KdMODIS_poids_IOCCG_NOMAD_BOUM_ss_12_COASTCOLOUR_443_to_670_log_Kd_lambda_merge_seuil_15_angle_ascii_6x1_hh_10_7_publi_CSIRO.sn" self.rinf_LUT_MOY = "Moy_KdMODIS_IOCCG_NOMAD_BOUM_ss_12_COASTCOLOUR_412_to_670_log_Kd_lambda_merge_seuil_15_ss_645_publi_inf_CSIRO.dat" self.rinf_LUT_ECART = "Ecart_KdMODIS_IOCCG_NOMAD_BOUM_ss_12_COASTCOLOUR_412_to_670_log_Kd_lambda_merge_seuil_15_ss_645_publi_inf_CSIRO.dat" # ======================================================================================================= # SEAWIFS if self.sensor == 'idSEAWIFS': self.WL = np.array([412, 443, 490, 510, 555, 670]) # Variables definiton for network Rrs490/Rrs555 >= .85 ------------------------------------------------- # hidden layer number self.rsup_NC=2 # input data number self.rsup_NE=5 # number of neuron in the first hidden layer self.rsup_NC1=4 # number of neuron in the second hidden layer self.rsup_NC2=4 # nb output self.rsup_NS=1 # nb input + nb output*/ self.rsup_NES=6 # LUTs file names self.rsup_LUT_POIDS = "KdSeaWiFS_poids_IOCCG_NOMAD_BOUM_ss_12_COASTCOLOUR_443_to_555_log_Kd_lambda_merge_seuil_15_angle_ascii_6x1_hh_4_4_publi_200715.sn" self.rsup_LUT_MOY = "Moy_KdSeaWiFS_IOCCG_NOMAD_BOUM_ss_12_COASTCOLOUR_412_to_670_log_Kd_lambda_merge_seuil_15_angle_publi_160715.dat" self.rsup_LUT_ECART = "Ecart_KdSeaWiFS_IOCCG_NOMAD_BOUM_ss_12_COASTCOLOUR_412_to_670_log_Kd_lambda_merge_seuil_15_angle_publi_160715.dat" # Variables definiton for network < .85 ----------------------------------- # nombre de couches cachees - hidden layer number self.rinf_NC=2 # input data number*/ self.rinf_NE=6 # number of neuron in the first hidden layer self.rinf_NC1=5 # number of neuron in the second hidden layer self.rinf_NC2=5 # nb output self.rinf_NS=1 # nb input + nb output self.rinf_NES=7 #LUTs file names self.rinf_LUT_POIDS = "KdSeaWiFS_poids_IOCCG_NOMAD_BOUM_ss_12_COASTCOLOUR_443_to_670_log_Kd_lambda_merge_seuil_15_angle_ascii_6x1_hh_5_5_publi_160715_ter.sn" self.rinf_LUT_MOY = "Moy_KdSeaWiFS_IOCCG_NOMAD_BOUM_ss_12_COASTCOLOUR_412_to_670_log_Kd_lambda_merge_seuil_15_angle_publi_160715.dat" self.rinf_LUT_ECART = "Ecart_KdSeaWiFS_IOCCG_NOMAD_BOUM_ss_12_COASTCOLOUR_412_to_670_log_Kd_lambda_merge_seuil_15_angle_publi_160715.dat" # LUTs arrays definition self.rsup_b1 = np.zeros((self.rsup_NC1), dtype=np.float64) self.rsup_b2 = np.zeros((self.rsup_NC2), dtype=np.float64) self.rsup_b3 = 0. self.rsup_w1 = np.zeros((self.rsup_NE,self.rsup_NC1), dtype=np.float64) self.rsup_w2 = np.zeros((self.rsup_NC1,self.rsup_NC2), dtype=np.float64) self.rsup_w3 = np.zeros((self.rsup_NC2), dtype=np.float64) self.rsup_moy = np.zeros((self.rsup_NES), dtype=np.float64) self.rsup_ecart = np.zeros((self.rsup_NES), dtype=np.float64) self.rinf_b1 = np.zeros((self.rinf_NC1), dtype=np.float64) self.rinf_b2 = np.zeros((self.rinf_NC2), dtype=np.float64) self.rinf_b3 = 0. self.rinf_w1 = np.zeros((self.rinf_NE,self.rinf_NC1), dtype=np.float64) self.rinf_w2 = np.zeros((self.rinf_NC1,self.rinf_NC2), dtype=np.float64) self.rinf_w3 = np.zeros((self.rinf_NC2), dtype=np.float64) self.rinf_moy = np.zeros((self.rinf_NES), dtype=np.float64) self.rinf_ecart = np.zeros((self.rinf_NES), dtype=np.float64) # -- METHOD read_LUTs_Kd ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ # Perceptron LUTs reading # LUTS are expected to be located in the directory defined by environment variable $KD_LUTS_PATH (cf: __init__) # -- def read_LUTs_Kd(self): print "class KdJamet.read_LUTs_Kd Method" # ---------------------------------------------------------------------------------- # < .85 weights settings #print " KdJamet.read_LUTs_Kd : read weights file for Rrs%d/Rrs%d < .85: %s/%s" % (self.WL[2], self.WL[4],self.lutpath, self.rinf_LUT_POIDS) weights = np.loadtxt( "%s/%s" % (self.lutpath, self.rinf_LUT_POIDS), skiprows=1) base = 0 self.rinf_b1 = weights[base:base+self.rinf_NC1,2] base = base+self.rinf_NC1 self.rinf_b2 = weights[base:base+self.rinf_NC2,2] base = base+self.rinf_NC2 self.rinf_b3 = weights[base:base+self.rinf_NS,2] base = base+self.rinf_NS for i in range(self.rinf_NE): self.rinf_w1[i,:] = weights[base:base+self.rinf_NC1,2] base = base+self.rinf_NC1 for i in range(self.rinf_NC1): self.rinf_w2[i,:] = weights[base:base+self.rinf_NC2,2] base = base+self.rinf_NC2 self.rinf_w3 = weights[base:base+self.rinf_NC2,2] usedData = np.array([1,2,3,4,5,6,8]) if self.sensor == 'idMODIS': usedData = np.array([1,2,3,4,6,7,9]) #print " KdJamet.read_LUTs_Kd : read mean file : %s/%s" % (self.lutpath, self.rinf_LUT_MOY) moy = np.loadtxt( "%s/%s" % (self.lutpath, self.rinf_LUT_MOY)) self.rinf_moy = moy[usedData] #print " class KdJamet.read_LUTs_Kd : read ecart file : %s/%s" % (self.lutpath, self.rinf_LUT_ECART) ecart = np.loadtxt( "%s/%s" % (self.lutpath, self.rinf_LUT_ECART)) self.rinf_ecart = ecart[usedData] # ---------------------------------------------------------------------------------- # >= .85 weights settings #print " KdJamet.read_LUTs_Kd : read weights file for Rrs%d/Rrs%d >= .85: %s/%s" % (self.WL[2], self.WL[4], self.lutpath, self.rsup_LUT_POIDS) weights = np.loadtxt( "%s/%s" % (self.lutpath, self.rsup_LUT_POIDS), skiprows=1) base = 0 self.rsup_b1 = weights[base:base+self.rsup_NC1,2] base = base+self.rsup_NC1 self.rsup_b2 = weights[base:base+self.rsup_NC2,2] base = base+self.rsup_NC2 self.rsup_b3 = weights[base:base+self.rsup_NS,2] base = base+self.rsup_NS for i in range(self.rsup_NE): self.rsup_w1[i,:] = weights[base:base+self.rsup_NC1,2] base = base+self.rsup_NC1 for i in range(self.rsup_NC1): self.rsup_w2[i,:] = weights[base:base+self.rsup_NC2,2] base = base+self.rsup_NC2 self.rsup_w3 = weights[base:base+self.rsup_NC2,2] usedData = np.array([1,2,3,4,6,8]) if self.sensor == 'idMODIS': usedData = np.array([1,2,3,4,7,9]) moy = np.loadtxt( "%s/%s" % (self.lutpath, self.rsup_LUT_MOY)) self.rsup_moy = moy[usedData] ecart = np.loadtxt( "%s/%s" % (self.lutpath, self.rsup_LUT_ECART)) self.rsup_ecart = ecart[usedData] # -- METHOD compute_allKd ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ # Kd computing at all Sensor WaveLength (vary for each sensor) # Inputs - # - indata : array[] [ Rrs443, Rrs490, Rrs510, Rrs555, Rrs665 ] # Output - # - return Kd[] at Sensor Wave Length # -- def compute_allKd(self, inRrs): outKd = np.zeros(NBWL,dtype=np.float32) if inRrs[2] / inRrs[4] < .85 : indata = np.zeros(self.rinf_NE,dtype=np.float32) indata[0:self.rinf_NE-1] = inRrs[1:self.rinf_NE] idLambda = self.rinf_NE-1 FLAGinf = True else: indata = np.zeros(self.rsup_NE,dtype=np.float32) indata[0:self.rsup_NE-1] = inRrs[1:self.rsup_NE] idLambda = self.rsup_NE-1 FLAGinf = False inc = 0 # compute Kd for each Wave Length of sensor used as iput for wl in self.WL: indata[idLambda] = wl if FLAGinf == True: outKd[inc] = self.passe_avant_inf(indata) else: outKd[inc] = self.passe_avant_sup(indata) inc = inc+1 return outKd # -- METHOD passe_avant_sup ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ # Kd computing for Rrs490/Rrs55 >= .85 at ONE WaveLength # Inputs - # - indata : array[rinf_NE] [ Rrs443, Rrs490, Rrs510, Rrs55, Lambda ] # Output - # - return y (Kd) at indata Wave Length {the wave length can be any value between min and max Wave length # used for the sensor } # -- def passe_avant_sup(self, indata): # Normalization x = ((2/3.) * (indata[:] - self.rsup_moy[0:self.rsup_NE])) / self.rsup_ecart[0:self.rsup_NE] a = np.zeros(self.rsup_NC1,dtype=np.float64) for i in range(self.rsup_NC1): for j in range(self.rsup_NE): a[i] = a[i] + x[j] * self.rsup_w1[j][i] a[:] = 1.715905*np.tanh((2./3.)*(a[:] + self.rsup_b1[:])) b = np.zeros(self.rsup_NC2,dtype=np.float64) for i in range(self.rsup_NC2): for j in range(self.rsup_NC1): b[i] = b[i] + a[j] * self.rsup_w2[j][i] b[:] = 1.715905*np.tanh((2./3.)*(b[:] + self.rsup_b2[:])) y = np.sum(b[:] * self.rsup_w3[:]) # Denormalisation y = 1.5*(y + self.rsup_b3)*self.rsup_ecart[self.rsup_NES-1] + self.rsup_moy[self.rsup_NES-1] y = 10.**y return(y) # -- METHOD passe_avant_inf ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ # Kd computing for Rrs490/Rrs55 < .85 at ONE WaveLength # Inputs # - indata : array[rinf_NE] [ Rrs443, Rrs490, Rrs510, Rrs55, Rrs665, Lambda ] # Output - # - return y (Kd) at indata Wave Length {the wave length can be any value between min and max Wave length # used for the sensor } # -- def passe_avant_inf(self, indata): # Normalization x = ((2/3.) * (indata[:] - self.rinf_moy[0:self.rinf_NE])) / self.rinf_ecart[0:self.rinf_NE] a = np.zeros(self.rinf_NC1,dtype=np.float64) for i in range(self.rinf_NC1): for j in range(self.rinf_NE): a[i] = a[i] + x[j] * self.rinf_w1[j][i] a[:] = 1.715905*np.tanh((2./3.)*(a[:] + self.rinf_b1[:])) b = np.zeros(self.rinf_NC2,dtype=np.float64) for i in range(self.rinf_NC2): for j in range(self.rinf_NC1): b[i] = b[i] + a[j] * self.rinf_w2[j][i] b[:] = 1.715905*np.tanh((2./3.)*(b[:] + self.rinf_b2[:])) y = np.sum(b[:] * self.rinf_w3[:]) # Denormalisation y = 1.5*(y + self.rinf_b3)*self.rinf_ecart[self.rinf_NES-1] + self.rinf_moy[self.rinf_NES-1] y = 10.**y return(y)