function DealiasedData = Dealiase(OriginalData)
% This script uses the following m-files:
    % alias_filter.m
    % spectral_functions.m

u = OriginalData; % The input data is re-named for simplicity
N = length(u);
dt = 3600;         % Because our data is hourly data = 3600 seconds (our preferred unit)
spectral_functions % Generates anonymous functions to calculate auto- or cross-spectra

% ---------------------------------------------------
%        Generation of the raw spectrum:
% ---------------------------------------------------
f    = fun.f(u,dt);
S    = fun.sxy(u,u,dt);
teta = angle(fft(fun.data(u)));

% ---------------------------------------------------
%   Correcting the spectra for the aliasing effect:
% ---------------------------------------------------
fc = 1;
f_limit = 1/(24*3600);
[temp.a,S_filtered,temp.b,temp.c] = alias_filter(f,S,1/dt,fc,f_limit,0);


% ---------------------------------------------------
%       Generating the de-aliased time series:
% ---------------------------------------------------
% This is done by reversing the de-aliased spectrum, keeping the phase of
% the original time series
X = sqrt(S_filtered*N/dt);       % since S = abs(fft(data)).^2*dt/N;
X = [0;X;X(end-1:-1:1)];         % since s = s(2:N/2+1) and s is symmetric
X = cos(teta).*X+i*sin(teta).*X; % the original phase

DealiasedData = real(ifft(X));           % The imaginary parts are negletably small
DealiasedData = DealiasedData + mean(u); % Adding the correct mean wind speed

