function out=hedgeMakeCloud2(data,xs)
%hedgeMakeCloud2 takes historical cloud data and generates a random
%sequence drawn from the same PDF. This is then shuffled until
%its autocorrelation is sufficiently similar to the autocorrelation of the
%original sequence.
%=========================================================================
% Usage:
% out=hedgeMakeCloud2(data,xs)
% Inputs:
% * data=vector of cloud cover fractions
% * xs = bins for the PDF, in NOAA cloud cover units (i.e. 0:10)
% Output, out = new random cloud cover vector in NOAA cloud cover units
%length of cloud cover vector
numDays=length(data);
%calculate target autocorrelation
xcorr1=xcorr(data,data);
% get PDF of values in cloud data
ys=hist(data,xs);
%populate a vector with random values from that PDF
dataTemp=zeros(numDays,1);
for i=1:numDays
dataTemp(i)=randarb(xs,ys);
end
% calculate autocorrelation of random vector as a baseline
score=xcorr(dataTemp,dataTemp);
%monte carlo search parameters:
scoreTemp=sqrt(sum((score-xcorr1).^2)); %current dist between autocorrelations
bestScore=scoreTemp;
closeEnough=1000; %what distance is good enough - this is empirical, for seasons of length 214 days.
%initialize tick counter for monte carlo
i=0;
%number of sequential days to move around at once
windowSize=2;
%shuffle until the score is good enough.
while scoreTemp>closeEnough
%find two random points in the vector to swap.
which1=ceil(rand()*(numDays-windowSize));
which2=ceil(rand()*(numDays-windowSize));
%switch the data (and subsequent window) at these points
dataProposed=dataTemp;
dataProposed(which1:which1+windowSize)=dataTemp(which2:which2+windowSize);
dataProposed(which2:which2+windowSize)=dataTemp(which1:which1+windowSize);
%update counter
i=i+1;
%get new distance score
xcorrTemp=xcorr(dataProposed,dataProposed);
scoreTemp=sqrt(sum((xcorrTemp-xcorr1).^2));
% if it's an improvement, keep it, and update monte carlo parameters
if scoreTemp < bestScore
bestScore=scoreTemp;
% disp(corr2(xcorrTemp,xcorr1)); %comment this line in to see correlation between autocorrelation vectors
dataTemp=dataProposed;
end
end
%output the final sequence
out=dataTemp;