diff utils/SegUtil.py @ 14:6dae41887406

added funcs for ssm enhancement in segutil scipt

--- a/utils/SegUtil.py	Fri Jun 05 18:02:05 2015 +0100
+++ b/utils/SegUtil.py	Fri Jun 12 17:44:39 2015 +0100
@@ -19,9 +19,11 @@
 from scipy.spatial import distance
 from scipy.ndimage import filters, zoom
 from scipy import signal
+from scipy.signal import correlate2d, convolve2d, filtfilt, resample, butter
 import pylab as plt
 from scipy.spatial.distance import squareform, pdist
-from scipy.ndimage.filters import *
+from scipy.ndimage.filters import maximum_filter, minimum_filter, percentile_filter, uniform_filter
+from scipy.ndimage.filters import median_filter as med_filter
 from sklearn.metrics.pairwise import pairwise_distances
 
 
@@ -45,7 +47,6 @@
 	if not os.path.exists(directory):
 		os.makedirs(directory)
 
-
 def median_filter(X, M=8):
 	"""Median filter along the first axis of the feature matrix X."""
 	for i in xrange(X.shape[1]):
@@ -293,6 +294,24 @@
 		return np.min(X)
 	fMin = 0
 	return fMin
+
+def lp(signal, fc=0.34, axis=-1):
+	'''Low pass filter function
+	signal: Raw signal to be smoothed.
+	fc: Cutoff frequency of the butterworth filter. Normalized from 0 to 1, where 1 is the Nyquist frequency.
+	axis: The axis of x to which the filter is applied. Default is -1.'''
+	bCoeffs, aCoeffs = butter(2, fc)
+	lp_smoothed_signal = filtfilt(bCoeffs, aCoeffs, signal, axis)
+	return lp_smoothed_signal
+
+def hp(signal, fc=0.34, axis=-1):
+	'''Low pass filter function
+	signal: Raw signal to be smoothed.
+	fc: Cutoff frequency of the butterworth filter.
+	axis: The axis of x to which the filter is applied. Default is -1.'''
+	bCoeffs, aCoeffs = butter(2, fc, 'highpass')
+	hp_smoothed_signal = filtfilt(bCoeffs, aCoeffs, signal, axis)
+	return hp_smoothed_signal
 	
 def getMean(feature, winlen, stepsize):
 	means = []
@@ -315,7 +334,7 @@
 	return delta_feature
 
 
-def getSSM(feature_array, metric='cosine', norm='simple', reduce=False):
+def getSSM(feature_array, metric='cosine', norm='exp', reduce=False):
 	'''Compute SSM given input feature array. 
 	args: norm: ['simple', 'remove_noise']
 	'''
@@ -323,13 +342,50 @@
 	dm = np.nan_to_num(dm)
 	if norm == 'simple':
 		ssm = 1 - (dm - np.min(dm)) / (np.max(dm) - np.min(dm))
+	if norm == 'exp': # Use with cosine metric only
+		ssm = np.exp(dm - 1)
 	if reduce:
 		ssm = reduceSSM(ssm)
 	return ssm
 
+def enhanceSSM(ssm, fc=0.34, med_size=(5,5), max_size=(5,5), min_size=(5,5), filter_type='min', axis=-1):
+	'''A series of filtering for SSM enhancement
+	fc: cutoff frequency for LP filtering.
+	med_size: Median filter window size.
+			  int or tuple. If using an integer for a 2d input, axis must be specified.
+	filter_type: Select either to use maximum filter or minimum filter.
+			float ['min', 'max', None]
+	max_size: Maximum filter window size.
+			  int or tuple. If using an integer for a 2d input, axis must be specified.
+			  Use this when homogeneity in the SSM is expressed by LARGE value.
+	min_size: Mininum filter window size.
+			  int or tuple. If using an integer for a 2d input, axis must be specified.
+			  Use this when homogeneity in the SSM is expressed by SMALL value. 
+			  (eg. When cosine metric and exp normalization and used for distance computation.)'''
 
+	ssm_lp = lp(enhanced_ssm, fc=fc)
+	
+	# Use scipy.ndimage.filters.median_filter instead
+	ssm_med = med_filter(ssm_lp, size=med_size)
+	
+	if filter_type == 'min':
+		enhanced_ssm = minimum_filter(ssm_med, size=min_size)
+	elif filter_type == 'max':
+		enhanced_ssm = maximum_filter(ssm_med, size=max_size)	
+	else:
+		enhanced_ssm = ssm_med
+	return enhanced_ssm
+	
 def reduceSSM(ssm, maxfilter_size = 2, remove_size=50):
-	reduced_ssm = ssm
+	'''Adaptive thresholding using OTSU method
+	Required package: skimage (0.10+)'''
+	
+	from skimage.morphology import disk
+	# from skimage.filters import threshold_otsu, rank #skimage 0.12
+	from skimage.filter.rank import otsu #skimage 0.10
+	from skimage.filter import threshold_otsu
+	
+	reduced_ssm = copy(ssm)
 	reduced_ssm[reduced_ssm<0.75] = 0
 	# # reduced_ssm = maximum_filter(reduced_ssm,size=maxfilter_size)
 	# # reduced_ssm = morphology.remove_small_objects(reduced_ssm.astype(bool), min_size=remove_size)