In [3]:
%pylab inline
#%matplotlib notebook
import glob
import os
import numpy as np
import scipy.signal as sig
from matplotlib import image as im
from matplotlib import pyplot as plt
from matplotlib.colors import LogNorm
import pywt 
#from pywt import thresholding as thr
#import nbsp


Populating the interactive namespace from numpy and matplotlib
---------------------------------------------------------------------------
ImportError                               Traceback (most recent call last)
<ipython-input-3-87fb6afd0851> in <module>()
      8 from matplotlib import pyplot as plt
      9 from matplotlib.colors import LogNorm
---> 10 import pywt
     11 #from pywt import thresholding as thr
     12 #import nbsp

ImportError: No module named 'pywt'

In [4]:
def load_PHD(path, query):
    files = glob.glob(path + query)
    files.sort()
    ascs = []
    filenames = []
    for i in files:
        with open(i, 'r') as f:
            FHD_size = asarray(f.readline().strip().split("\t"),dtype=int)    
        filenames.append(os.path.basename(i))                  
        ascs.append(swapaxes(swapaxes(np.genfromtxt(i, skip_header=1).reshape(FHD_size[2],FHD_size[0],FHD_size[1]),0,2),0,1))
    images_arr = np.asarray(ascs)
    print("Archivos cargados: " + str(images_arr.shape[0]))
    print(filenames)
    return images_arr, filenames

In [5]:
from scipy import ndimage
from numpy import fft

def fft_gauss (in_array, kernel_size, array_size):
    
    # FFT filter
    im_fft = fft.rfftn(in_array)
    im_rfft_filtered = ndimage.fourier_gaussian(im_fft, kernel_size, array_size[1])
    im_filtered = fft.irfftn(im_rfft_filtered)
    
    #Power spectrum       
    pwr_spectrum = abs(fft.fftshift(im_fft))**2
    pwr_spectrum_filtered = abs(fft.fftshift(im_rfft_filtered))**2
        
    #Re-normalization
    #for i in np.arange(in_array.shape[0]):
    sum_ratio= in_array.sum()/im_filtered.sum()
    im_filtered = im_filtered*sum_ratio
    
    return im_filtered, pwr_spectrum, pwr_spectrum_filtered

#Vector version
vfft_gauss = np.vectorize(fft_gauss, excluded = ['kernel_size', 'array_size'], otypes = [np.ndarray] )

In [25]:
path = "/home/nicolas/Documentos/Universidad/Simulaciones/Jnirs/newCode/IM/"
PHD_fluo, PHD_filenames = load_PHD(
    path, "*300.dat_PHD-Ascii.dat")

print(PHD_fluo[0].shape)


Archivos cargados: 1
['jnirs-IM-z1s0_a300.dat_PHD-Ascii.dat']
(100, 100, 25)

In [27]:
plt.imshow(PHD_fluo[0][:,:,10])
cbar = plt.colorbar()
#FHD_fondo_fluo[FHD_fondo_filenames.index("fondo-fast2_FHD-Ascii.dat")][:,:,0]



In [ ]:


In [42]:
div1=FHD_fluo[FHD_filenames.index("inclusion-vicky-c_FHD-Ascii.dat")]/FHD_fondo_fluo[FHD_fondo_filenames.index(fondo)]
div1_fft=fft_gauss(div1[40,:,:], 1, div1.shape)
plt.imshow(div1_fft[0])
cbar = plt.colorbar()



In [9]:
div1=FHD_fluo[FHD_filenames.index("inclusion-x21-z115_FHD-Ascii.dat")]/FHD_fondo_fluo[FHD_fondo_filenames.index("fondo-fl_FHD-Ascii.dat")]
div1_fft=fft_gauss(div1[:,:,6], 1, div1.shape)
plt.imshow(div1_fft[0], vmax=1.2, vmin=0.9)
cbar = plt.colorbar()


/usr/local/lib/python2.7/dist-packages/ipykernel/__main__.py:1: RuntimeWarning: divide by zero encountered in divide
  if __name__ == '__main__':

In [18]:
divf=FHD_fondo_fluo[FHD_fondo_filenames.index("fondo-slow-fl_FHD-Ascii.dat")]/FHD_fondo_fluo[FHD_fondo_filenames.index("fondo-fl_FHD-Ascii.dat")]
plt.imshow(divf[50,:,:], vmax=1.1)
cbar = plt.colorbar()


/usr/local/lib/python2.7/dist-packages/ipykernel/__main__.py:1: RuntimeWarning: divide by zero encountered in divide
  if __name__ == '__main__':

In [ ]: