function [spiketrain, isi] = inhPoisson(DT,r,f,T)

     % Function that returns a point process generated by the density lambda_t
     % IN: timestep DT, rate r, frequency f, Duration T

m=0;  % time index
tmp_spiketrain(1) = -Inf;   % temporary variable for the spiketrain
n=1;

% Loop

while  n==1 | tmp_spiketrain(n-1) < T,
  x=-log(rand);
  I=0;
  mlocal=0;  % local index
  while I<x,
    I = I+DT*lambda_t(r,f,m*DT);
    m = m+1;
    mlocal = mlocal + 1;
  end
  
  if n>1
    tmp_isi(n-1) =  mlocal*DT;
    tmp_spiketrain(n) = tmp_spiketrain(n-1) + tmp_isi(n-1);
  else
    tmp_spiketrain(n) = mlocal*DT;
  end
  n=n+1;
end

% delte the last spike beyond the period T
spiketrain = tmp_spiketrain(1:end-1); 
isi = tmp_isi(1:end-1);