Initial commit.

checkModel.py

@@ -0,0 +1,69 @@
+import os
+import sys
+
+import h5py
+import numpy as np
+
+from scipy.io import savemat
+from dataLoaders import loadInputs
+from gammanet import createArchitecture
+
+n = int( sys.argv[1] )
+filePath = sys.argv[2]
+mode = sys.argv[3]
+
+if mode == "cars":
+    modelDir = "./models/cars/"
+else:
+    modelDir = "./models/raman/"
+
+dirs = os.listdir( modelDir )
+dirs = sorted( dirs )
+nModels = len( dirs )
+
+nInputs = 640
+nOutputs = 640
+
+modelPath = modelDir + dirs[nModels - 1]
+print( modelPath )
+
+model = createArchitecture( nInputs, nOutputs)
+
+file = h5py.File( modelPath, 'r')
+weights = []
+nWeights = len( file.keys() )
+
+for ii in range( nWeights ):
+   weights.append( file['weight' + str(ii)][:] )
+
+model.set_weights( weights )
+model.summary()
+
+X = loadInputs( filePath, mode)
+inputDims = (1, nInputs, 1)
+
+nDataSets = X.shape[0]
+result = {}
+
+nDataSets = min( n, nDataSets)
+
+for ii in range(nDataSets):
+
+    input = np.empty( inputDims )
+    input[ 0, :, 0] = X[ ii, :, 0]
+
+    yHat = model( input )
+
+    median = yHat.quantile( 0.50 )
+    lowerBound = yHat.quantile( 0.05 )
+    upperBound = yHat.quantile( 0.95 )
+
+    tempIndex = "spectrum_" + str(ii)
+
+    result[ tempIndex ] = {}
+    result[ tempIndex ]["input"] = input
+    result[ tempIndex ]["median"] = median.numpy()
+    result[ tempIndex ]["lowerBound"] = lowerBound.numpy()
+    result[ tempIndex ]["upperBound"] = upperBound.numpy()
+
+    savemat("./modelCheck.mat", result)

dataLoaders.py

@@ -0,0 +1,84 @@
+import h5py
+import numpy as np
+
+def loadData( fileName, mode, maximumDataSets = 1000000):
+
+    modes = {"raman", "cars"}
+
+    if mode not in modes:
+        raise ValueError('Mode must be either "raman" or "cars". Current: %s.', mode)
+
+    data = h5py.File( fileName, 'r')
+
+    measurementData = data.get( mode + "Data" )
+    measurementData = np.array( measurementData )
+
+    imChi3Data = data.get('beta')
+    imChi3Data = np.array( imChi3Data )
+
+    nPoints = measurementData.shape[0]
+    nDataSets = measurementData.shape[1]
+
+    if( nDataSets > maximumDataSets ):
+        nDataSets = maximumDataSets
+
+    xDims = ( nDataSets, nPoints, 1)
+    yDims = ( nDataSets, nPoints)
+
+    X = np.empty( xDims )
+    y = np.empty( yDims )
+
+    for ii in range( nDataSets ):
+
+        X[ ii, :, 0] = measurementData[ :, ii]
+        y[ ii, :] = imChi3Data[ :, ii]
+
+    return X, y
+
+def loadInputs( fileName, mode, maximumDataSets = 1000000):
+
+    modes = {"raman", "cars"}
+
+    if mode not in modes:
+        raise ValueError('Mode must be either "raman" or "cars". Current: %s.', mode)
+
+    data = h5py.File( fileName, 'r')
+
+    measurementData = data.get( mode + "Data" )
+    measurementData = np.array( measurementData )
+
+    nPoints = measurementData.shape[0]
+    nDataSets = measurementData.shape[1]
+
+    if( nDataSets > maximumDataSets ):
+        nDataSets = maximumDataSets
+
+    xDims = ( nDataSets, nPoints, 1)
+
+    X = np.empty( xDims )
+
+    for ii in range( nDataSets ):
+        X[ ii, :, 0] = measurementData[ :, ii]
+
+    return X
+
+def loadDataField( fileName, field, maximumDataSets = 1000000):
+
+    data = h5py.File( fileName, 'r')
+
+    fieldData = data.get( field )
+    fieldData = np.array( fieldData )
+
+    nOutputPoints = fieldData.shape[0]
+    nDataSets = fieldData.shape[1]
+
+    if( nDataSets > maximumDataSets ):
+        nDataSets = maximumDataSets
+
+    yDims = ( nDataSets, nOutputPoints)
+    y = np.empty( yDims )
+
+    for ii in range( nDataSets ):
+            y[ ii, :] = fieldData[ :, ii]
+
+    return y

gammanet.py

@@ -0,0 +1,58 @@
+import keras
+from keras.layers import Activation, BatchNormalization, Conv1D, Dense, Dropout, Flatten
+from keras.models import Sequential
+
+import tensorflow as tf
+import tensorflow_probability as tfp
+
+tfd = tfp.distributions
+tfpl = tfp.layers
+
+optimizer = keras.optimizers.Adam( learning_rate = 0.001 )
+
+def createArchitecture( inputSize, outputSize):
+
+    inputShape = ( inputSize, 1)
+
+    negativeLogLikelihood = lambda y, rv_y: -rv_y.log_prob( y )
+    gammaLayer = lambda input: tfd.Gamma( concentration = tf.nn.relu( input[..., 0:outputSize] + 1e-12 ),
+                                          log_rate      = tf.nn.relu( input[..., outputSize:]) )
+
+    outputLayer = tfp.layers.DistributionLambda( gammaLayer )
+
+    model = Sequential()
+    model.add( BatchNormalization(  axis = -1,
+                                    momentum = 0.99,
+                                    epsilon = 0.001,
+                                    center = True,
+                                    scale = True,
+                                    beta_initializer = 'zeros',
+                                    gamma_initializer = 'ones',
+                                    moving_mean_initializer = 'zeros',
+                                    moving_variance_initializer = 'ones',
+                                    beta_regularizer = None,
+                                    gamma_regularizer = None,
+                                    beta_constraint = None,
+                                    gamma_constraint = None,
+                                    input_shape = inputShape
+                                 ))
+
+    model.add( Activation('relu'))
+
+    model.add( tfp.layers.Convolution1DFlipout( 128, activation = 'relu', kernel_size = (32)) )
+    model.add( Conv1D( 64, activation = 'relu', kernel_size = (16)) )
+    model.add( Conv1D( 16, activation = 'relu', kernel_size = (8)) )
+    model.add( Conv1D( 16, activation = 'relu', kernel_size = (8)) )
+    model.add( Conv1D( 16, activation = 'relu', kernel_size = (8)) )
+
+    model.add( Dense(  32, activation = 'relu') )
+    model.add( Dense(  16, activation = 'relu') )
+    model.add( Flatten() )
+    model.add( Dropout( 0.25 ) )
+
+    model.add( Dense( outputSize + outputSize, activation = 'relu' ) )
+    model.add( outputLayer )
+
+    model.compile( loss = negativeLogLikelihood, optimizer = optimizer)
+
+    return model

trainNetwork.py

@@ -0,0 +1,75 @@
+import sys
+import h5py
+import numpy as np
+
+from datetime import datetime
+from scipy.io import savemat
+from sklearn.model_selection import train_test_split
+
+from gammanet import createArchitecture
+from dataLoaders import loadData
+
+import pandas as pd
+import tensorflow as tf
+
+tf.keras.backend.clear_session()
+tf.compat.v1.enable_eager_execution()
+
+fileName = sys.argv[1]
+mode = sys.argv[2]
+epochs = int( sys.argv[3] )
+batch_size = int( sys.argv[4] )
+validationSplit = float( sys.argv[5] )
+
+timeStamp = datetime.now()
+timeStamp = datetime.timestamp( timeStamp )
+timeStamp = int( timeStamp )
+timeStamp = str( timeStamp )
+
+def scheduler( epoch, lr):
+    if epoch < 10: return lr
+    elif lr > 0.000001: return lr * tf.math.exp(-0.1)
+    else: return 0.000001
+
+
+callbackSchedule = tf.keras.callbacks.LearningRateScheduler(scheduler)
+callbackStopping = tf.keras.callbacks.EarlyStopping( monitor = 'val_loss',
+                                             patience = 50, 
+                                             min_delta = 0, 
+                                             mode = 'auto', 
+                                             baseline = None, 
+                                             restore_best_weights = True)
+
+X, y = loadData( fileName, mode)
+X_train, X_test, y_train, y_test = train_test_split( X, y, test_size = validationSplit)
+
+nInputs = X_train.shape[1]
+nOutputs = y_train.shape[1]
+
+model_ii = createArchitecture( nInputs, nOutputs)
+model_ii.summary()
+
+h = model_ii.fit( X_train, y_train, epochs = epochs,
+                    verbose = 1,
+                    validation_data = ( X_test, y_test),
+                    batch_size = batch_size,
+                    shuffle = True,
+                    callbacks = [callbackSchedule, callbackStopping])
+
+modelSavePath = './models/' + mode + '/gamma-specnet-model-' + timeStamp
+historySavePath = './models/' + mode + '/gamma-specnet-history-' + timeStamp + ".json"
+
+hist_df = pd.DataFrame( h.history ) 
+
+# save to json:  
+with open( historySavePath, mode = 'w') as f:
+    hist_df.to_json( f )
+
+file = h5py.File( modelSavePath, 'w')
+
+weights = model_ii.get_weights()
+nWeights = len( weights )
+
+for ii in range( nWeights ):
+   file.create_dataset( 'weight' + str(ii), data = weights[ii])
+file.close()