function [AP_vector, presynaptic_stimulus_matrix_tally]=random_train_AP_times_only(pulse_duration,samples_per_ms,samples_per_pulse,trial,...
    presynaptic_stimulus_matrix,presynaptic_stimulus_matrix_tally,sum_of_maximum_IPIs,postsynaptic_Vm_step,OVSF,ext_command_sensitivity)
%makes the train of ones and zeros to evoke APs in the presynaptic cell
%except for 'cannot occur' related lines, proofread nov 2 2006 -- appears to require about 3ms per iteration

row_one=presynaptic_stimulus_matrix(1,:);   %the first row of presynaptic_stimulus_matrix is the number of possible values for each IPI

%cannot_occur=[3 5 6 7 8 10 11 12 13 14 15 16 17 18]   %THIS LINE is part of the effort to reduce the possible values of IPIa
if trial==1
    %na=cannot_occur(1)                     %THIS LINE is part of the effort to reduce the possible values of IPIa
    %while find (na(1)== cannot_occur)      %THIS LOOP is part of the effort to reduce the possible values of IPIa 
        na=rand(1,length(row_one));
        na=na.*row_one;
        na=ceil(na);                        %then na goes and provides the info needed to create the output train
    %end
    presynaptic_stimulus_matrix_tally=zeros(row_one);
    if length(row_one)<2                    %the eventuality where row_one is only one element long, that is, there is only one IPI, that is, there are only two presynaptic APs
        presynaptic_stimulus_matrix_tally=zeros(row_one,1);
    end
    dummest=length(row_one);
    nummest=prod(row_one);
    while dummest                           %this loop exists so that this emmfile can handle presynaptic_stimulus_matrix of any size (that is, any number of IPIs)
        nummest=nummest-(((row_one(dummest)-na(dummest))/row_one(dummest))*prod(row_one(1:dummest)));
        dummest=dummest-1;    
    end
    presynaptic_stimulus_matrix_tally(nummest)=presynaptic_stimulus_matrix_tally(nummest)+1;
    %presynaptic_stimulus_matrix_tally(cannot_occur)=1                  %THIS LINE is part of the effort to reduce the possible values of IPIa
else
    min_before=min(presynaptic_stimulus_matrix_tally);
    for i=1:length(size(presynaptic_stimulus_matrix_tally))-1         %this tomfoolery is necessary bc 'min' takes the min along each dimension, so need 'number of [non-singleton] dimensions of presynaptic_stimulus_matrix_tally' times
        min_before=min(min_before);           %the first time is outside the for loop (if one or more of presynaptic_stimulus_matrix_tally's dimensions are singleton, this won't change the loop's behavior)  
    end
    pseudomin=min_before-1;                 %this is just done just to start the while loop on the right track
    while min_before~=pseudomin
        na=rand(1,length(row_one));
        na=na.*row_one;
        na=ceil(na); 
        dummest=length(row_one);
        nummest=prod(row_one);
        while dummest                       %this loop exists so that this emmfile can handle presynaptic_stimulus_matrix of any size (that is, any number of IPIs)
            nummest=nummest-(((row_one(dummest)-na(dummest))/row_one(dummest))*prod(row_one(1:dummest)));
            dummest=dummest-1;    
        end
        pseudomin=presynaptic_stimulus_matrix_tally(nummest);            %escape the while loop if the randomly generated value's spot in the presynaptic_stimulus_matrix_tally matrix has the same value as the min of the presynaptic_stimulus_matrix_tally matrix
    end
    presynaptic_stimulus_matrix_tally(nummest)=presynaptic_stimulus_matrix_tally(nummest)+1;
    %if min(presynaptic_stimulus_matrix_tally(:,1))==max(presynaptic_stimulus_matrix_tally(:,1))    %THIS LOOP was added for when each possible IPI has happened exactly as many times as each other possible IPI -- part of the "THIS LINE" modification oct 29 2003
        %presynaptic_stimulus_matrix_tally(cannot_occur)=presynaptic_stimulus_matrix_tally(cannot_occur)+1   %THIS LINE is part of the effort to reduce the possible values of IPIa
    %end
end

ba=na-1;                                    %correcting, for indexing (see next few lines)
intervuls=ba.*presynaptic_stimulus_matrix(3,:);
intervuls=intervuls+presynaptic_stimulus_matrix(2,:);
intervuls=intervuls*samples_per_ms;         %translating the intervals from milliseconds into samples
intervuls=[0 intervuls];
spiketimes=cumsum(intervuls);

AP_vector=zeros(sum_of_maximum_IPIs,1);
%each pulse lasts from the assigned time + one sample to the assigned time + a pulse's worth of samples
%example: if the pulse is set for 600 ms, sampling 10 kHz, and pulse_duration 1 ms, the pulse will last from samples 6001-6010, inclusive
for j=1:length(spiketimes)                     %the length of ba = the number of IPIs -- with the initial pulse @ halfsecond, the # of IPIs = the # of remaining pulses
    AP_time=spiketimes(j);                     %now interval is beginning of one pulse to beginning of the next (old way was end of one to beginning of the next -- 
    AP_vector(AP_time+1:AP_time+samples_per_pulse)=1;
    %AP_vector(AP_time+1:AP_time+floor(samples_per_pulse/2))=-(postsynaptic_Vm_step/OVSF);               %if the *presyn* Vm is modulated, Na activation can be addressed by down-up pulses
    %AP_vector(AP_time+1+ceil(samples_per_pulse/2):AP_time+samples_per_pulse)=((40-(-70+(postsynaptic_Vm_step*ext_command_sensitivity)))/ext_command_sensitivity)/OVSF;    %ideal is to have the 'down' part back to -70, the 'up' up to ~ +30
end