Biological Machines/Sensory Systems/Visual System/Simulation of the Visual System

From Wikibooks, open books for an open world
Jump to: navigation, search

In this section an overview in the simulation of processing done by the early levels of the visual system will be given. The implementation to reproduce the action of the visual system will thereby be done with MATLAB and its toolboxes. The processing done by the early visual system was discussed in the section before and can be put together with some of the functions they perform in the following schematic overview. A good description of the image processing can be found in (Cormack 2000).

Schematic overview of the processing done by the early visual system
Structure Operations 2D Fourier Plane
World I(x,y,t,\lambda) 2D Fourier Plane 01.jpg
Optics Low-pass spatial filtering 2D Fourier Plane 02.jpg
Photoreceptor Array Sampling, more low-pass filtering, temporal lowhandpass filtering, \lambda filtering, gain control, response compression
LGN Cells Spatiotemporal bandpass filtering, \lambda filtering, multiple parallel representations 2D Fourier Plane 03.jpg
Primary Visual Cortical Neurons: Simple & Complex Simple cells: orientation, phase, motion, binocular disparity, & \lambda filtering 2D Fourier Plane 04.jpg
Complex cells: no phase filtering (contrast energy detection)

On the left, are some of the major structures to be discussed; in the middle, are some of the major operations done at the associated structure; in the right, are the 2-D Fourier representations of the world, retinal image, and sensitivities typical of a ganglion and cortical cell. (From Handbook of Image and Video Processing, A. Bovik)

As we can see in the above overview different stages of the image processing have to be considered to simulate the response of the visual system to a stimulus. The next section will therefore give a brief discussion in Image Processing. But first of all we will be concerned with the Simulation of Sensory Organ Components.

Simulating Sensory Organ Components[edit]

Anatomical Parameters of the Eye[edit]

The average eye has an anterior corneal radius of curvature of r_C = 7.8 mm , and an aqueous refractive index of 1.336. The length of the eye is L_E = 24.2 mm. The iris is approximately flat, and the edge of the iris (also called limbus) has a radius r_L = 5.86 mm.

Optics of the Eyeball[edit]

The optics of the eyeball are characterized by its 2-D spatial impulse response function, the Point Spread Function (PSF)


h(r) = 0.95\cdot \exp\left( -2.6\cdot |r|^{1.36} \right) + 0.05\cdot\exp\left( -2.4\cdot |r|^{1.74} \right) ,

in which r is the radial distance in minutes of arc from the center of the image.

Practical implementation[edit]

Obviously, the effect on a given digital image depends on the distance of that image from your eyes. As a simple place-holder, substitute this filter with a Gaussian filter with height 30, and with a standard deviation of 1.5.

In one dimension, a Gaussian is described by


g(x) = a \cdot \exp \left( -\frac{x^2}{2\sigma^2} \right) .

Activity of Ganglion Cells[edit]

Ignoring the

  • temporal response
  • effect of wavelength (especially for the cones)
  • opening of the iris
  • sampling and distribution of photo receptors
  • bleaching of the photo-pigment

we can approximate the response of ganglion cells with a Difference of Gaussians (DOG, Wikipedia [1])


f(x;\sigma) = \frac{1}{\sigma_1\sqrt{2\pi}} \, \exp \left( -\frac{x^2}{2\sigma_1^2} \right)-\frac{1}{\sigma_2\sqrt{2\pi}} \, \exp \left( -\frac{x^2}{2\sigma_2^2} \right).

The values of \sigma_1 and \sigma_2 have a ratio of approximately 1:1.6, but vary as a function of eccentricity. For midget cells (or P-cells), the Receptive Field Size (RFS) is approximately

RFS \approx 2 \cdot \text{Eccentricity} ,

where the RFS is given in arcmin, and the Eccentricity in mm distance from the center of the fovea (Cormack 2000).

Activity of simple cells in the primary visual cortex (V1)[edit]

Again ignoring temporal properties, the activity of simple cells in the primary visual cortex (V1) can be modeled with the use of Gabor filters (Wikipedia [2]). A Gabor filter is a linear filter whose impulse response is defined by a harmonic function (sinusoid) multiplied by a Gaussian function. The Gaussian function causes the amplitude of the harmonic function to diminish away from the origin, but near the origin, the properties of the harmonic function dominate

g(x,y;\lambda,\theta,\psi,\sigma,\gamma)=\exp\left(-\frac{x'^2+\gamma^2y'^2}{2\sigma^2}\right)\cos\left(2\pi\frac{x'}{\lambda}+\psi\right) ,

where

x' = x \cos\theta + y \sin\theta\, ,

and

y' = -x \sin\theta + y \cos\theta\, .

In this equation, \lambda represents the wavelength of the cosine factor, \theta represents the orientation of the normal to the parallel stripes of a Gabor function (Wikipedia [3]), \psi is the phase offset, \sigma is the sigma of the Gaussian envelope and \gamma is the spatial aspect ratio, and specifies the ellipticity of the support of the Gabor function.

Gabor function, with sigma = 1, theta = 1, lambda = 4, psi = 2, gamma = 1

This is an example implementation in MATLAB:

function gb = gabor_fn(sigma,theta,lambda,psi,gamma)
 
  sigma_x = sigma;
  sigma_y = sigma/gamma;
 
  % Bounding box
  nstds = 3;
  xmax = max(abs(nstds*sigma_x*cos(theta)),abs(nstds*sigma_y*sin(theta)));
  xmax = ceil(max(1,xmax));
  ymax = max(abs(nstds*sigma_x*sin(theta)),abs(nstds*sigma_y*cos(theta)));
  ymax = ceil(max(1,ymax));
  xmin = -xmax;
  ymin = -ymax;
  [x,y] = meshgrid(xmin:xmax,ymin:ymax);
 
  % Rotation
  x_theta = x*cos(theta) + y*sin(theta);
  y_theta = -x*sin(theta) + y*cos(theta);
 
  gb = exp(-.5*(x_theta.^2/sigma_x^2+y_theta.^2/sigma_y^2)).* cos(2*pi/lambda*x_theta+psi);
 
end