function PAR=PARabalone4(Nvar,Nspecies,Nx,Ny)
%
%   define model parameters and values
%    some parameters are independent of space, others are arrays
PAR.case='4';

PAR.Nvar=Nvar;
PAR.Nspecies=Nspecies;
PAR.Nx=Nx;
PAR.Ny=Ny;
PAR.iS=1;PAR.iI=2;PAR.iD=3;PAR.iP=4;

PAR.depth = ones(Ny,Nx);   %   vertical dimension of the volume (m)

% infection rate by contact with infections particles
%  [infected produced/infectious particle/day]
%  the convention is IPinfect(source species, infected species)
PAR.IPinfect = zeros(Nspecies,Nspecies);  
%  particles from species 1 infect species 2
PAR.IPinfect(1,2) = 0.025;
PAR.IPinfect(2,1) = 0.001;
%PAR.IPinfect(1,1) = .003;
%PAR.IPinfect(2,2) = 0.025;
%PAR.IPinfect = 0.025*ones(Nspecies,Nspecies);  

% infection rate by contact with infectious individuals 
%  [infected produced/susceptible animal/day]
PAR.Iinfect = zeros(Nspecies,Nspecies);  
%PAR.Iinfect = 0.02*ones(Nspecies,Nspecies);  

% infection rate by contact with dead infected individuals 
%  [infected produced/susceptible animal/day]
PAR.Dinfect = zeros(Nspecies,Nspecies);  
%PAR.Dinfect = 0.0008*ones(Nspecies,Nspecies);  

%  reproduction rate for susceptibles
PAR.Srepro = zeros(Nspecies,1);  %
PAR.Srepro(1,:,:) = .01;
PAR.Srepro(2,:,:) = .01;
%  reproduction rate for infected
PAR.Irepro = zeros(Nspecies,1);  
%PAR.Irepro(1) = .05;
%PAR.Irepro(2) = .05;
%   carrying capacity for the population
PAR.Carry = zeros(Nspecies,1);  
PAR.Carry(1)=150;
PAR.Carry(2)=100;

% mortality rate of infected [1/day]
PAR.Imort = 0.08*ones(Nspecies,1);  

% background mortality rate of infected [1/day]
%PAR.Bmort = 0.08*ones(Nspecies,1);  
PAR.Bmort = zeros(Nspecies,1);  

% removal rate of dead [1/day]
PAR.DeadDecay = 1.5*ones(Nspecies,1);  

% infectious particles released by infected [particles/animal/day]
%PAR.Irelease  = .015*ones(Nspecies,1); 
PAR.Irelease  = zeros(Nspecies,1); 
PAR.Irelease(1)=.1;
PAR.Irelease(2)=.01;
% infectious particles released by dead [particles/animal/day]
PAR.Drelease = ones(Nspecies,1);  

% removal rate of IP from the environment [1/day]
%PAR.IPremove = 0.001*ones(Nspecies,1);  
PAR.IPremove = zeros(Nspecies,1);  
PAR.IPremove(1) = 0.001;
PAR.IPremove(2) = 0.001;

% rate of transfer of IP in the E/W direction [1/day]
%PAR.Uex = zeros(Ny,Nx-1);  
%PAR.Uex = 0.05*ones(Ny,Nx-1);  
%PAR.Uex = [0.02 0.02 ;0 0; -0.02 -0.02 ];  
PAR.Uex = [0.02 0.02 ;0 0;0 0; -0.02 -0.02 ];  

% rate of transfer of IP in the N/S direction [1/day]
%PAR.Vex = zeros(Ny-1,Nx);  
%PAR.Vex = 0.2*ones(Ny-1,Nx);  
%PAR.Vex = [-0.02 0  0.02; -0.02 0  0.02];
PAR.Vex = [-0.02 0  0.02; -0.02 0  0.02;  -0.02 0  0.02];