function spikes = lifboltzmannspikes(trials, input, tmax, gain)
% Generate spike times of a leaky integrate-and-fire neuron.
% trials: the number of trials to be generated.
% input: the stimulus intensity.
% tmax: the duration of a trial.
% gain: gain of the neuron, i.e. the slope factor of the boltzmann input.

    tau = 0.01;
    D = 1e-2;
    imax = 25;
    ithresh = 10;
    slope = gain;
    vreset = 0.0;
    vthresh = 10.0;
    dt = 1e-4;
    
    n = ceil(tmax/dt);
    inb = imax./(1.0+exp(-slope.*(input - ithresh)));
    spikes = cell(trials, 1);
    for k=1:trials
        times = [];
        j = 1;
        v = vreset + (vthresh - vreset) * rand();
        noise = sqrt(2.0*D)*randn(n, 1)/sqrt(dt);
        for i=1:n
            v = v + (- v + noise(i) + inb)*dt/tau;
            if v >= vthresh
                v = vreset;
                times(j) = i*dt;
                j = j + 1;
            end
        end
        spikes{k} = times;
    end
end