rand_int.m

我认为很不错的语音处理的ｍａｔｌａｂ源代码

        n=prod(n);
    else
        mm=n;
        nn=1;
    end

    % Make sure b > a
    % swap values if necessary
    if a > b
        c=a;
        a=b;
        b=c;
    end

    % Make sure b-a >= 1 and that n >= 1.
    % Make sure there are at least n integers from b to a.
    % There are b-a+1 integers from b to a inclusive.
    % if conditions are false return reasonable output.
    if logical((b-a+1) < 1)
        flag=1;
        ndraw=[];
        warning('Input error!  Empty range of integers.  Returning empty output array');
    end

    % If specified to select less than 1 integer return an empty matrix.
    if logical(n < 1)
        flag=1;
        ndraw=[];
        warning('Input error!  Select less than 1 integer.  Returning empty output array');
    end

    % If specified to select more integers than available then return all
    % availabe integers and print a warning.
    if logical(n > (b-a+1)) && isequal(no_duplicates, 1)
        flag=1;
        ndraw=(a:b)';
        warning('Input error!  Selecting more integers than in the range and choosing no duplicates.  Returning a:b as the output array');
    end

    % if flag is set to 1 then set the error
    % and truncate the output to reasonable output.
    if isequal(flag, 1)
        error=1;
        n2=min([length(ndraw), n]);
        if n2 >= 1
            ndraw=ndraw(1:n2, 1);
        else
            ndraw=[];
        end
    end

    n3=n;
    if n >= (b-a+1)/2
        %reverse algorithm and deselect data points
        flag2=1;
        n=(b-a+1)-n;
    end

    % define the threshold limits for the different regimes
    thresh_regime1=0.43;    % generates entire matrix

    % The regime 2 threshold depends on whether selecting almost all the
    % data points or very few data points
    % not sure what teh optimum threshold is
    % this does not crash
    thresh_regime2l=min([0.07, 7/n3]);  % The probability of a duplicate
    % grows with range and number
    % of selected integers
    thresh_regime2h=min([0.05, 5/((b-a+1)-n3)]);  %

    if min((b-a+1)-n3, n3)/(b-a) > thresh_regime1
        % Regime 1
        % n is nearly equal to half the number of data points
        flag3=1;
    elseif logical(n3 <= (b-a+1)*thresh_regime2l) || logical(n3 >= (b-a+1)*(1-thresh_regime2h))
        % Regime 2
        % n is a small fraction of the range
        flag3=2;
    else
        % Regime 3 catch any failures from the first two regimes
        % Regime 2 will fail. regime 1 should b never fail.
        % Select random integers then remove the extra integers
        flag3=3;
    end

    % maximum number of times the selection while loop will run
    max_count=max(ceil(n3/(n3-(b-a))), 100);

    % nunique is the number of unique integers selected
    nunique=0;

    % count is the number of iterations of the selection while loop
    count=0;

    % count2 is the number of iterations of the deselection while loop
    count2=0;

    num_points=b-a+1;

    if isequal(flag, 1)

    else
        if isequal(no_duplicates, 1)

            % For the case of removing all duplicates
            fail=0;

            if  isequal(flag3, 1)
                % In this regime the entire matrix of values if filled, 
                % then the extras are discarded.
                % This is a faster approach for the case
                % min((b-a+1)-n,n)/(b-a) > thresh_regime1
                % Code was contributed by Jos x@y.z
                % code was modified by sorting the output
                % when sorting is specified
                r=randperm(num_points);
                if isequal(sorted, 1)
                    ndraw=sort(a-1+r(1:n))';
                else
                    ndraw=a-1+r(1:n)';
                end

            elseif isequal(flag3, 2)
                % This is the 100 percent purge regime.
                % There is a possibility of returning an empty array.
                %
                % The program runs faster by rejecting all integers rather than
                % sorting through the input and keeping the good integers
                %
                % This is faster approach for the case
                % (b-a+1)/n > 1-thresh_regime2 || n/(b-a) < thresh_regime2
                %
                % The judicious cuttoff of ten percent minimizes the chances of
                % an error and returning an empty array.

                flag5 = 1;
                count=0;
                while flag5 && logical(count < max_count)
                    count=count+1;
                    ndraw = unique(round((a - 1/2) + (b-a+1)*rand(n,1)))';
                    flag5 = (length(ndraw)~=n);
                end

            end

            % check to find out if enough integers were selected
            if logical(length(ndraw)~=n) && ~isequal(flag3, 3)
                ndraw=[];
                error=1;
                fail=1;
                % This warning is caught by the if statement below
                % warning('Failed to select enough integers');
            end

            % if not enough integers were selected
            if (~isequal(flag3, 1) && ~isequal(flag3, 2)) || fail

                ndraw=[];
                count=0;

                while (logical(nunique < n) || logical(count < 1)) && logical(count < max_count)

                    count=count+1;

                    % To make all combinations more equally likely, stretch the random
                    % distribution to a-0.5 and b+0.5.
                    % Round the numbers to the nearest integers.
                    % Make sure the selected integers are unique.
                    % Only one selection per integer.
                    % select an extra thirty percent of points because there will
                    % duplicates that are discarded.
                    num_draws=max(ceil((1.1+n/num_points)*(n-nunique)));
                    ndraw1=round((a+b)/2+((b+0.5)-(a-0.5))*(-0.5+rand(num_draws, 1)));

                    % make sure the numbers fit in the proper range
                    ndraw1(ndraw1 > b)=b;
                    ndraw1(ndraw1 < a)=a;

                    ndraw=unique([ndraw' ndraw1'])';

                    % Count the number of unique integers selected
                    nunique=length(ndraw);
                    count2=0;

                    if nunique > n
                        ndraw_del=[];
                        num_draws=ceil(1.3*(nunique-n));
                        aa=1;
                        bb=nunique;
                        nunique_del=0;

                        while nunique_del < (nunique-n) && logical(count2 < 100)
                            count2=count2+1;
                            ndraw_del=unique([ndraw_del' round((aa+bb)/2+((bb+0.5)-(aa-0.5))*(-0.5+rand(num_draws, 1)))']');
                            nunique_del=length(ndraw_del);
                        end
                        p = randperm(nunique_del);
                        delete_ix=ndraw_del(p);

                        save_ix=setdiff(aa:bb, delete_ix(1:(nunique-n)));
                        ndraw=ndraw(save_ix);

                    end

                end

            end


            % If the number of integers to be selected is greater than half the
            % size of the draw pool of then use a deselection method.
            if isequal(flag2, 1)
                if length(ndraw) < 2
                    ndraw=setdiff((a:b)', ndraw');
                else
                    ndraw=setdiff((a:b)', ndraw')';
                end
            end

            if ~isequal(sorted, 1)
                draw_pool=randperm(length(ndraw))';
                ndraw=ndraw(draw_pool);
            end

            if fail
                if length(ndraw)~=n3
                    warning('Returning empty array');
                    ndraw=[];
                    error=1;
                else
                    % The error has been succesfully corrected 
                    %warning('Error Canceled.  Correct Number of Integers Selected.  Returned full array');
                    error=0;
                end
            end

        else

            % For the case of allowing duplicates

            if isequal(flag2, 1)
                n=(b-a+1)-n;
            end

            ndraw=round((a+b)/2+((b+0.5)-(a-0.5))*(-0.5+rand(n, 1)));

            % make sure the numbers fit in the proper range
            ndraw(ndraw > b)=b;
            ndraw(ndraw < a)=a;

            if isequal(sorted, 1)
                % sort output
                ndraw=sort(ndraw);
            end

        end

    end
end

%
% ************************************************************************
%
% Reshape the output array to the desired dimensions
% If there are not enough data points, then
% determine the best size for the output.
%
if logical(nn >= 1) && logical(mm >= 1)
    if logical( numel(ndraw) >= mm*nn)
        % reshape the data to the specified size
        ndraw=reshape(ndraw, mm, nn);
    else
        % Not enough data points to
        % reshape the data to the specified size.
        %
        % Get the number of data points available
        % and reshape the output to a reasonable size.
        np=numel(ndraw);
        if np < 1
            ndraw=[];
        else

            % Find the smaller dimension size
            [md min_dim]=min([mm, nn]);

            % determine the best way to reshape the output matrix
            % try to preserve the value of the smaller dimension
            if isequal(min_dim, 1)
                if mm <= np
                    nn=floor(np/mm);
                else
                    mm=np;
                    nn=1;
                end
            else

                if nn <= np
                    mm=floor(np/nn);
                else
                    nn=np;
                    mm=1;
                end
            end

            ndraw=reshape(ndraw(1:(mm*nn)), mm, nn);
            
        end

    end

end

% Plot the output if desired
% Good for inspecting the randomness of the output.
if isequal(make_plot, 1) && ~isempty(ndraw)
    
    figure(1);
    delete(1);
    figure(1);
    [m1, n1]=size(ndraw);
    
    if m1 > n1
        ndraw2=ndraw';
        axis_dir='Row';
    else
        ndraw2=ndraw;
        axis_dir='Column';
    end
    
    [m1, n1]=size(ndraw2);
    
    plot(1:n1, ndraw2, 'linestyle', 'none', 'Marker', '.', 'MarkerEdgeColor', 'k');
    ylim([a b]);
    
    xlim([0.5 n1+0.5]);
    
    regime_desc={'Build Entire Array, Then Discard Extras', 'Pure Purge', 'Select 10% Extra ,Then Discard Extras'};
    
    title({'Plot of Random Integers', ['Regime ', num2str(flag3), ' ', regime_desc{flag3}]}, 'Fontsize', 18);
    xlabel([axis_dir,' Indices of Integers'], 'Fontsize', 16);
    ylabel('Value of Integers', 'Fontsize', 16);
    
end