function spikes = lifspikes(trials, input, tmaxdt, D)
% Generate spike times of a leaky integrate-and-fire neuron.
% trials: the number of trials to be generated
% input: the stimulus either as a single value or as a vector
% tmaxdt: in case of a single value stimulus the duration of a trial
%         in case of a vector as a stimulus the time step
% D: the strength of additive white noise
    
    tau = 0.01;
    if nargin < 4
        D = 1e0;
    end
    vreset = 0.0;
    vthresh = 10.0;
    dt = 5e-5;
    
    if max(size(input)) == 1
        input = input * ones(ceil(tmaxdt/dt), 1);
    else
        dt = tmaxdt;
    end
    spikes = cell(trials, 1);
    for k=1:trials
        times = [];
        j = 1;
        v = vreset + (vthresh-vreset)*rand(1);
        noise = sqrt(2.0*D)*randn(length(input), 1)/sqrt(dt);
        for i=1:length(noise)
            v = v + (- v + noise(i) + input(i))*dt/tau;
            if v >= vthresh
                v = vreset;
                spiketime = i*dt;
                times(j) = spiketime;
                j = j + 1;
            end
        end
        spikes{k} = times;
    end
end