%!PS-Adobe-2.0 EPSF-2.0
%%Title: 2-D Gabor function
%%Creator: Izumi Ohzawa, Ph.D. izumi@pinoko.berkeley.edu
%%CreationDate: Thu Oct  2 15:21:07 PDT 1997
%%Pages: 0 0
%%BoundingBox: 0 0 216 216
%%Orientation: Portrait
%%EndComments

% ==========================================================================================
% Copyright 1997, Izumi Ohzawa.
% You may freely distribute and modify this PostScript program for any
% purpose as long as this copyright notice and the following references
% are retained in the distribution.

% Note that Y and X dimensions are parallel and perpendicular, respectively,
% to grating bars.  Since grating bars are horizontal, Y is horizontal,
% and X is vertical.  Rotate as you see fit when importing the EPS.
% 

% ==========================================================================================
% Configurable Parameters
% ==========================================================================================

% ##### Parameter definitions
/CT_peak 1.0 def	% contrast at center
/Ncyc	 4.0 def	% spatial freq in # of cycles per domain (fractional values OK)
/Phase   0.0 def	% phase angle in degrees (0=Even sym., 90.0=Odd sym.)
/Ky	 2.3 def	% Inverse of Gaussian sigma for length dimension
/Kx	 2.3 def	% Inverse of Gaussian sigma for width dimension

% ##### Number of steps in X and Y dimensions
/NX 108 def
/NY 108 def

% Standard length unit definitions
/inch { 72.0 mul } bind def		% for obstinate Americans
/cm { 28.3464567 mul } bind def		% centimeter for all others

% ##### Paper/Display/Film dimension definitions
%  If you change these, change %%BoundingBox: at the top.
/Xsize { 3.0 inch } bind def			% width
/Ysize { 3.0  inch } bind def			% height

% ==========================================================================================
%                        No user serviceable part below this line.
% ==========================================================================================

/E	2.718281828 def
/Xstripwidth Xsize NY div def		% length of grating bars
/NX_2 NX 2 div def
/NY_2 NY 2 div def
/DataString NX string def

%%EndProlog
 
% %%%%%% Execution begins here %%%%%%

0 1 NY 1 sub {
    NY_2 sub			% index - (NY/2)
    NY_2 div			% -1 .. +1
    Ky mul dup mul neg E exch exp
    /CT_y exch def			% Y Envelope

    /currentSF Ncyc NX div  def	% Spatial frequency

    % Fill DataString with swept-frequency sine wave data
    0 1 NX 1 sub {
	dup
	NX_2 sub
	NX_2 div					% X = -1 .. +1
	Kx mul dup mul neg E exch exp			% exp(- K * X^2)
	CT_peak mul CT_y mul /CT exch def

	DataString exch	dup		% string index index
	NX_2 sub currentSF mul	% string index ((index-NX/2)*currentSF)
	360 mul Phase add cos		% string index cos((index-NX/2)*currentSF*360+Phase)
	127.4 mul			% string index 127*cos(...)
	CT mul 128 add cvi		% string index 128+ Contrast*127*cos(...)
	put				% store value in DataString at index
    } for

    gsave
    Xstripwidth Ysize scale		% map unit image to this size rect
    1 NX 8 [1 0 0 NX 0 0] { DataString } image	% this lays down bitmap
    grestore

    Xstripwidth 0 translate		% move over to next strip position
} for

showpage         			% is allowed in EPSF
%%Trailer