From c897257e5ab27c1a6ab104e6579bf02ebe5fa352 Mon Sep 17 00:00:00 2001 From: Tbs_Fchnr Date: Fri, 30 Oct 2020 10:35:19 +0100 Subject: [PATCH] Add fourier and histogram filter classes. Add notebook for histogram exploration. Implemented histogram equalisation functionality. --- TF00-initial-exploration.ipynb | 2 +- TF01-fourier-exploration.ipynb | 186 +++++++++++++++++++++++++++++++++ filters.py | 154 +++++++++++++++++++++++++-- handler.py | 26 +++-- requirements.txt | 3 + 5 files changed, 353 insertions(+), 18 deletions(-) create mode 100644 TF01-fourier-exploration.ipynb diff --git a/TF00-initial-exploration.ipynb b/TF00-initial-exploration.ipynb index 0643d36..b7818fd 100644 --- a/TF00-initial-exploration.ipynb +++ b/TF00-initial-exploration.ipynb @@ -106,7 +106,7 @@ "metadata": {}, "outputs": [], "source": [ - "gaus = filters.Gaussian(sig=2)\n", + "gaus = filters.Gaussian(sig=3)\n", "print(gaus)" ] }, diff --git a/TF01-fourier-exploration.ipynb b/TF01-fourier-exploration.ipynb new file mode 100644 index 0000000..baaee28 --- /dev/null +++ b/TF01-fourier-exploration.ipynb @@ -0,0 +1,186 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Fourier Exploration" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Imports" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "import filters\n", + "import handler" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "import numpy as np # for debugging purposes" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "import modules for debugging" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "from math import ceil, floor" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Load Image" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "img, _ = handler.getImageAsArray(handler.FOETUS_PATH_ORIGINAL)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "fourierFilter = filters.FFT_TruncateCoefficients(keep=0.1)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "imgFFT = fourierFilter.compute(img, plot=True)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "handler.plotFigs([img, imgFFT])" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "imgFFT2 = np.arange(25).reshape((5,5))\n", + "imgFFT2" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# Get shape of image: rows and columns\n", + "row, col = imgFFT2.shape\n", + "print(row, col)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "keep = 0.1\n", + "print(keep)\n", + "print(1-keep)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "print('imgFFT2[{}:{}, :]'.format(ceil(row*keep), floor(row*(1-keep))))\n", + "print('imgFFT2[:, {}:{}]'.format(ceil(col*keep), floor(col*(1-keep))))" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# Set all rows and cols to zero not within the keep fraction\n", + "imgFFT2[ceil(row*keep):floor(row*(1-keep)), :] = 0\n", + "imgFFT2[:, ceil(col*keep):floor(col*(1-keep))] = 0" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "imgFFT2" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.8.1" + } + }, + "nbformat": 4, + "nbformat_minor": 4 +} diff --git a/filters.py b/filters.py index 4aa557e..58c94ca 100644 --- a/filters.py +++ b/filters.py @@ -7,10 +7,13 @@ from abc import ABC, abstractmethod import numpy as np import matplotlib.pyplot as plt -from math import ceil +from math import ceil, floor import statistics +from scipy import fftpack +from matplotlib.colors import LogNorm -class Filter(ABC): +class SpatialFilter(ABC): + # TODO: implement standard mask shapes of square or cross and implement kernel creation based on this for each filter def __init__(self, maskSize, kernel, name, linearity): self.assertTypes(maskSize, kernel) self.name = name @@ -80,7 +83,82 @@ def convolve(self, img, padding=True): return output -class Median(Filter): +class FourierFilter: + def fft2D_scipy(self, img, plot=False): + """ + Function transforms image into Fourier domain + :param plot: bool to configure plotting of fourier spectum. default=False + :param img: image to be transformed + :return: image in fourier domain/ fourier spectrum of image + """ + imgFFT = fftpack.fft2(img) + if plot: self.plotFourierSpectrum(imgFFT) + return imgFFT + + @staticmethod + def dft(x): # not in use + """ + Function computes the discrete fourier transform of a 1D array + :param x: input array, 1 dimensional + :return: np.array of fourier transformed input array + """ + x = np.asarray(x, dtype=float) + N = x.shape[0] + n = np.arange(N) + k = n.reshape((N, 1)) + M = np.exp((-2j * np.pi * k * n) / N) + return np.dot(M, x) + + def fft(self, x): + """ + Function recursively implements 1D Cooley-Turkey fast fourier transform + :param x: input array, 1 dimensional + :return: np.array of fourier transformed input array + """ + x = np.array(x, dtype=float) + N = x.shape[0] + + if N % 2 > 0: + raise ValueError("size of x must be a power of 2") + elif N <= 32: + return self.dft(x) + else: + X_even = self.fft(x[::2]) + X_odd = self.fft(x[1::2]) + factor = np.exp((-2j * np.pi * np.arange(N)) /N) + return np.concatenate([X_even + factor[:N / 2] * X_odd, + X_even + factor[N / 2:] * X_odd ]) + + + def fft2D(self, x): + """ + Function recursively implements 1D Cooley-Turkey fast fourier transform + :param x: input array, 1 dimensional + :return: np.array of fourier transformed input array + """ + x = np.array(x, dtype=float) + xRot = x.T + + self.fft(x) + + @staticmethod + def inverseFFT_scipy(img): + return fftpack.ifft2(img).real + + @staticmethod + def plotFourierSpectrum(imgFFT): + """ + Function displays fourier spectrum of image that has been fourier transformed + :param imgFFT: fourier spectrum of img + :return: None + """ + plt.figure() + plt.imshow(np.abs(imgFFT), norm=LogNorm(vmin=5)) + plt.colorbar() + plt.title('Fourier Spectrum') + + +class Median(SpatialFilter): def __init__(self, maskSize): kernel = np.zeros((maskSize,maskSize)) @@ -92,7 +170,7 @@ def __init__(self, maskSize): def computePixel(self, sub): return statistics.median(sub.flatten()) -class Mean(Filter): +class Mean(SpatialFilter): """ Effectively a low pass filter. Alternative kernel implemented in class LowPass(Filter). """ @@ -108,7 +186,7 @@ def computePixel(self, sub): # element-wise multiplication of the kernel and image pixel under consideration return (self.kernel * sub).sum() -class Gaussian(Filter): +class Gaussian(SpatialFilter): def __init__(self, sig): # Calculate mask size from sigma value. Ensures filter approaches zero at edges (always round up) maskSize = ceil((6 * sig) + 1) @@ -124,7 +202,7 @@ def __init__(self, sig): kernel = np.exp(-0.5 * (np.square(xx) + np.square(yy)) / np.square(sig)) - super().__init__(maskSize, kernel, name='gaussian', linearity='non-linear') + super().__init__(maskSize, kernel, name='gaussian', linearity='linear') def computePixel(self, sub): """ @@ -135,7 +213,10 @@ def computePixel(self, sub): """ return (self.kernel * sub).sum()/ self.kernel.sum() -class HighPass(Filter): +class HighPass(SpatialFilter): + """ + High pass filter to have sharpening effect on image. + """ def __init__(self, maskSize): # TODO: Make ratio of intensity reduction vs. increase configurable for both high and low pass @@ -153,7 +234,7 @@ def computePixel(self, sub): return (self.kernel * sub).sum() -class LowPass(Filter): +class LowPass(SpatialFilter): def __init__(self, maskSize): kernel = np.zeros((maskSize, maskSize)) @@ -165,4 +246,59 @@ def __init__(self, maskSize): super().__init__(maskSize, kernel, name='low-pass', linearity='non-linear') def computePixel(self, sub): - return (self.kernel * sub).sum()/ self.kernel.sum() \ No newline at end of file + return (self.kernel * sub).sum()/ self.kernel.sum() + +class FFT_TruncateCoefficients(FourierFilter): + def __init__(self, keep=0.1): + self.keep = keep + + def compute(self, img, plot=False): + # Get fourier transform of image + imgFFT = self.fft2D_scipy(img, plot=plot) + + # Call ff a copy of original transform + imgFFT2 = imgFFT.copy() + + # Get shape of image: rows and columns + row, col = imgFFT2.shape + + # Set all rows and cols to zero not within the keep fraction + imgFFT2[ceil(row*self.keep):floor(row*(1-self.keep)), :] = 0 + imgFFT2[:, ceil(col*self.keep):floor(col*(1-self.keep))] = 0 + + if plot: self.plotFourierSpectrum(imgFFT2) + + return self.inverseFFT_scipy(imgFFT2) + +class HistogramFilter: + @staticmethod + def getHistogram(img): + # Create zeros array with as many items as bins to collect for + histogram = np.zeros(256) + + # Loop through image and add pixel intensities to histogram + for pixel in img.flatten(): + histogram[pixel] += 1 + + cumSum = HistogramFilter.cumSum(histogram) + cumSumNormalised = HistogramFilter.cumSumNormalise(cumSum) + + return histogram, cumSumNormalised + + @staticmethod + def cumSum(histogram): + return np.cumsum(histogram.flatten()) + + @staticmethod + def cumSumNormalise(cumSum): + # normalize the cumsum to range 0-255 + csNormalised = 255 * (cumSum - cumSum.min()) / (cumSum.max() - cumSum.min()) + + # cast it to uint8 + return csNormalised.astype('uint8') + + @staticmethod + def histEqualise(img, cs): + imgNew = cs[img.flatten()] + + return np.reshape(imgNew, img.shape) \ No newline at end of file diff --git a/handler.py b/handler.py index 69bca11..82a4536 100644 --- a/handler.py +++ b/handler.py @@ -18,8 +18,6 @@ FOETUS_PATH_ORIGINAL = ".\\images\\foetus.png" NZJERS_PATH_ORIGINAL = ".\\images\\NZjers1.png" -FOETUS_PATH_FILTERED = ".\\images\\foetus-filtered.png" -NZJERS_PATH_FILTERED = ".\\images\\NZjers1-filtered.png" def getImageAsArray(path, convertToGrey=True): """ @@ -58,17 +56,28 @@ def plotFigs(images): else: raise Exception("Make sure you pass in either a single image as np ndarray or list of images as np ndarray.") - for img in images: - # Convert array to geyscale pillow image object - img_PIL = Image.fromarray(img, 'L') - display(img_PIL) + # set up side-by-side image display + fig = plt.figure() + fig.set_figheight(15) + fig.set_figwidth(15) + + ax1 = fig.add_subplot(1,2,1) + ax1.title.set_text('original') + plt.imshow(images[0], cmap='gray') + + # display the new image + ax2 = fig.add_subplot(1,2,2) + ax2.title.set_text('filtered') + plt.imshow(images[1], cmap='gray') + + plt.show(block=True) def saveAll(img, filtr, saveFilter=True): """ Function to save all figures relevant to report. Currently filtered image and plot of kernel. :return: """ - assert isinstance(filtr, filters.Filter) + assert isinstance(filtr, filters.SpatialFilter) currentDir = Path().absolute() root = str(currentDir) + '\\..\outputs\{}\maskSize_{}\\'.format(filtr.name, filtr.maskSize) @@ -99,4 +108,5 @@ def saveAll(img, filtr, saveFilter=True): f.write(''.join("filter.{} = {}\n".format(k, v))) print("Saved filter object attributes to... \n{}\n\n".format(root + 'filter.txt')) else: - pass \ No newline at end of file + pass + diff --git a/requirements.txt b/requirements.txt index 3533b51..377c9eb 100644 --- a/requirements.txt +++ b/requirements.txt @@ -16,6 +16,7 @@ ipython-genutils==0.2.0 ipywidgets==7.5.1 jedi==0.17.2 Jinja2==2.11.2 +joblib==0.17.0 jsonschema==3.2.0 jupyter==1.0.0 jupyter-client==6.1.7 @@ -49,11 +50,13 @@ pywinpty==0.5.7 pyzmq==19.0.2 qtconsole==4.7.7 QtPy==1.9.0 +scikit-learn==0.23.2 scipy==1.5.2 Send2Trash==1.5.0 six==1.15.0 terminado==0.9.1 testpath==0.4.4 +threadpoolctl==2.1.0 tornado==6.0.4 traitlets==5.0.5 wcwidth==0.2.5