function freq_grid = gen_freqgrid(N, Fs, opt)
% gen_freqgrid   Generate frequency grid
%   freq_grid = gen_freqgrid(N, Fs, opt) generate an N point
%   frequency grid according to sample rate Fs. This grid matches
%   the operation used in fftshift(opt = 1) or not(opt = 0).
%
%   % Example:
%   %   Create a 16 point frequency grid for a sample rate of 10 
%   %   Hz used in fftshift.
%   
%   freq_grid = gen_freqgrid(16, 10, 1)
    
% set 'CenterDC' in psdfreqvec to true preserves Nyquist point,
% which does not match our processing to the data since we use
% fftshift.

% adopted from psdfreqvec 

%% Checking 'opt' in input
    if ~(opt == 0 || opt == 1)
        disp('Option must be 0(No fftshift) or 1(fftshift)');
        return;
    end

%% Generating freqeuncy grid
    % freq_grid = fftshift(psdfreqvec(...
    %     'Npts',N,'Fs',Fs,'Range','whole'));
    freq_res = Fs/N;
    freq_grid = (0:N-1).'*freq_res;

    if opt == 1
        % linspace(freq_offset-fs/2, freq_offset+fs/2*(fft_len-2)/fft_len, fft_len);
        Nyq = Fs/2;
        half_res = freq_res/2;
        if rem(N,2) % odd
            idx = 1:(N-1)/2;
            halfpts = (N+1)/2;
            freq_grid(halfpts) = Nyq-half_res;
            freq_grid(halfpts+1) = Nyq+half_res;
        else
            idx = 1:N/2;
            hafpts = N/2+1;
            freq_grid(hafpts) = Nyq;
        end
        freq_grid(N) = Fs-freq_res;
        freq_grid = fftshift(freq_grid);
        freq_grid(idx) = freq_grid(idx)-Fs;        
    end    
end