diff --git a/code/local_features.py b/code/local_features.py deleted file mode 100755 index 11b5b3b8f8e1ecdb0cd2941d0dc5b050d040a4d6..0000000000000000000000000000000000000000 --- a/code/local_features.py +++ /dev/null @@ -1,117 +0,0 @@ -""" -Helping functions for extracting features from images -""" - -import skimage.io -import numpy as np -import matplotlib.pyplot as plt -import scipy.ndimage - -def get_gauss_feat_im(im, sigma=1, normalize=False): - """Gauss derivative feaures for every image pixel. - Arguments: - image: a 2D image, shape (r,c). - sigma: standard deviation for Gaussian derivatives. - normalize: flag indicating normalization of features. - Returns: - imfeat: a 3D array of size (r,c,15) with a 15-dimentional feature - vector for every image pixel. - Author: vand@dtu.dk, 2020 - """ - - r,c = im.shape - imfeat = np.zeros((r,c,15)) - imfeat[:,:,0] = scipy.ndimage.gaussian_filter(im,sigma,order=0) - imfeat[:,:,1] = scipy.ndimage.gaussian_filter(im,sigma,order=[0,1]) - imfeat[:,:,2] = scipy.ndimage.gaussian_filter(im,sigma,order=[1,0]) - imfeat[:,:,3] = scipy.ndimage.gaussian_filter(im,sigma,order=[0,2]) - imfeat[:,:,4] = scipy.ndimage.gaussian_filter(im,sigma,order=[1,1]) - imfeat[:,:,5] = scipy.ndimage.gaussian_filter(im,sigma,order=[2,0]) - imfeat[:,:,6] = scipy.ndimage.gaussian_filter(im,sigma,order=[0,3]) - imfeat[:,:,7] = scipy.ndimage.gaussian_filter(im,sigma,order=[1,2]) - imfeat[:,:,8] = scipy.ndimage.gaussian_filter(im,sigma,order=[2,1]) - imfeat[:,:,9] = scipy.ndimage.gaussian_filter(im,sigma,order=[3,0]) - imfeat[:,:,10] = scipy.ndimage.gaussian_filter(im,sigma,order=[0,4]) - imfeat[:,:,11] = scipy.ndimage.gaussian_filter(im,sigma,order=[1,3]) - imfeat[:,:,12] = scipy.ndimage.gaussian_filter(im,sigma,order=[2,2]) - imfeat[:,:,13] = scipy.ndimage.gaussian_filter(im,sigma,order=[3,1]) - imfeat[:,:,14] = scipy.ndimage.gaussian_filter(im,sigma,order=[4,0]) - - if normalize: - imfeat -= np.mean(imfeat, axis=(0,1)) - imfeat *= (1/np.std(imfeat, axis=(0,1))) - - return imfeat - -def im2col(im, patch_size=[3,3], stepsize=1): - """Rearrange image patches into columns - Arguments: - image: a 2D image, shape (r,c). - patch size: size of extracted paches. - stepsize: patch step size. - Returns: - patches: a 2D array which in every column has a patch associated - with one image pixel. For stepsize 1, number of returned column - is (r-patch_size[0]+1)*(c-patch_size[0]+1) due to bounary. The - length of columns is pathc_size[0]*patch_size[1]. - """ - - r,c = im.shape - s0, s1 = im.strides - nrows =r-patch_size[0]+1 - ncols = c-patch_size[1]+1 - shp = patch_size[0],patch_size[1],nrows,ncols - strd = s0,s1,s0,s1 - - out_view = np.lib.stride_tricks.as_strided(im, shape=shp, strides=strd) - return out_view.reshape(patch_size[0]*patch_size[1],-1)[:,::stepsize] - - -def ndim2col(im, block_size=[3,3], stepsize=1): - """Rearrange image blocks into columns for N-D image (e.g. RGB image)""""" - if(im.ndim == 2): - return im2col(im, block_size, stepsize) - else: - r,c,l = im.shape - patches = np.zeros((l*block_size[0]*block_size[1], - (r-block_size[0]+1)*(c-block_size[1]+1))) - for i in range(l): - patches[i*block_size[0]*block_size[1]:(i+1)*block_size[0]*block_size[1], - :] = im2col(im[:,:,i],block_size,stepsize) - return patches - -#%% - -if __name__ == '__main__': - - filename = '../example_data/training_image.png' - I = skimage.io.imread(filename) - I = I.astype(np.float) - - - # features based on gaussian derivatives - sigma = 1; - gf = get_gauss_feat_im(I, sigma) - - fig,ax = plt.subplots(3,5) - for j in range(5): - for i in range(3): - ax[i][j].imshow(gf[:,:,5*i+j]) - ax[i][j].set_title(f'layer {5*i+j}') - - - # features based on image patches - I = I[I.shape[0]//2:I.shape[0]//2+20,I.shape[1]//2:I.shape[1]//2+20] # smaller image such that we can see the difference in layers - pf = im2col(I,[3,3]) - pf = pf.reshape((9,I.shape[0]-2,I.shape[1]-2)) - - fig, ax = plt.subplots() - ax.imshow(I) - - fig,ax = plt.subplots(3,3) - for j in range(3): - for i in range(3): - ax[i][j].imshow(pf[3*i+j]) - ax[i][j].set_title(f'layer {3*i+j}') - - diff --git a/code/microscopy_analysis.py b/code/microscopy_analysis.py index 394ea5f76edd06f7306ab44c7451b725f6e4a5f8..f5505e708cf7af40176f2aa9b1815c3f4f555cf4 100755 --- a/code/microscopy_analysis.py +++ b/code/microscopy_analysis.py @@ -10,7 +10,6 @@ import os import skimage.io as io import skimage.transform import numpy as np -import local_features as lf import matplotlib.pyplot as plt from math import ceil @@ -310,37 +309,6 @@ def rgb2cmy(img): #%%Functions for feature analysis and visualisation -# computes the max projection image and the features into a list -def compute_max_img_feat(in_dir, ext, base_name, sigma, sc_fac, save = False, abs_intensity = True): - dir_list = [dI for dI in os.listdir(in_dir) if (os.path.isdir(os.path.join(in_dir,dI)) and dI[0:len(base_name)]==base_name)] - dir_list.sort() - max_img_list = [] - for d in dir_list: - image_dir = in_dir + d + '/' - max_img = get_max_img(image_dir, ext) - - if save: - dir_name_out_img = in_dir.replace('data','maxProjImages') - os.makedirs(dir_name_out_img, exist_ok = True) - - io.imsave( dir_name_out_img + d + '.png', max_img.astype('uint8')) - - max_img_list += [skimage.transform.rescale(max_img, sc_fac, multichannel=True)] - - feat_list = [] - for max_img in max_img_list: - r,c = max_img.shape[:2] - feat = np.zeros((r*c,45)) - for i in range(0,3): - feat_tmp = lf.get_gauss_feat_im(max_img[:,:,i], sigma) - feat[:,i*15:(i+1)*15] = feat_tmp.reshape(-1,feat_tmp.shape[2]) - if not(abs_intensity): - feat_list += [np.concatenate((feat[:,1:15],feat[:,16:30],feat[:,31:45]),axis = 1)] - else: - feat_list += [feat] - return max_img_list, feat_list - - # computes a histogram of features from a kmeans object def compute_assignment_hist(feat_list, kmeans, background_feat = None): assignment_list = []