function [x] = generate_vector(N, K, ensemble, sigma, p)
% =========================================================================
%                       Sparse vector generator
% =========================================================================
% INPUT ARGUMENTS:
% N                         Size of sparse vector.
% K                         Sparsity of vector.
% ensemble                  Ensemble type of measurement matrix. Possible
%                           values are:
%                           -'Gaussian': K non-zero elements of the sparse vector
%                           are drawn from normal distribution N(0,1).
%                           -'sGaussian': K non-zero elements of the sparse vector
%                           are drawn from normal distribution N(0,sigma^2).
%                           -'Bernoulli': K non-zero elements of the sparse vector
%                           are drawn from Bernoulli distribution (1/2,1/2).
%                           -'pBernoulli': K non-zero elements of the sparse vector
%                           are drawn from Bernoulli distribution (p,1-p).
% sigma                     Standard deviation of Gaussian distribution.
% p                         Parameter of Bernoulli distribution.
% =========================================================================
% OUTPUT ARGUMENTS:
% x                     K-sparse vector.
% =========================================================================
% 01/04/2011, by Anastasios Kyrillidis. anastasios.kyrillidis@epfl.ch, EPFL.
% =========================================================================

if nargin < 3
    ensemble = 'Gaussian';
end;

if nargin < 4
    sigma = 1;
    p = 0.5; 
end;

x = zeros(N,1);
rand_indices = randperm(N);

switch ensemble
    case 'Gaussian'
        x(rand_indices(1:K)) = randn(K,1);         % Standard normal distribution ~ N(0,1)
    case 'sGaussian'
        x(rand_indices(1:K)) = sigma*randn(K,1);    % Normal distribution ~ N(0,sigma^2)
    case 'Uniform'
        x(rand_indices(1:K)) = rand(K,1);           % Uniform [0,1] distribution
    case 'Bernoulli'
        x(rand_indices(1:K)) = (-1).^round(rand(K,1));     % Bernoulli ~ (1/2, 1/2) distribution
    case 'pBernoulli'
        x(rand_indices(1:K)) = (-1).^(rand(K,1) > p);   % Bernoulli ~ (p, 1-p) distribution     
end;

x = x/norm(x);