%this script takes daily data for two timeseries of ensemble data as input and computes lagged correlations up to lags of 365 days
%the input timeseries are assumed to have dimension (nT,nEns)
function [lagged_corr,lagged_corr_detrend,lagged_corr_init,lagged_corr_detrend_init,lags]=compute_lagged_corr_daily_ensemble(predictand,predictor,varname1,varname2)

if(mod(size(predictand,1),365)~=0 || mod(size(predictor,1),365)~=0 )
disp('Error: Input timeseries are not an integer number of years')
return
end

nEns=size(predictand,2);

nTarget=365;
nLag=365;

lags=0:1:nLag-1;

lagged_corr=zeros(nTarget,nLag);
lagged_corr_detrend=zeros(nTarget,nLag);

shift=nLag;

for day=1:nTarget

day;

start_ind=shift+day;
target_tmp=predictand(start_ind:365:end,:);

for lag=0:nLag-1
predictor_tmp=predictor(start_ind-lag:365:end-lag,:);


r=corrcoef(target_tmp(:),predictor_tmp(:));
lagged_corr(day,lag+1)=r(1,2);

target_tmp_detrend=detrend(target_tmp);
predictor_tmp_detrend=detrend(predictor_tmp);

r=corrcoef(target_tmp_detrend(:),predictor_tmp_detrend(:));
lagged_corr_detrend(day,lag+1)=r(1,2);

end

end

nInit=365;
lagged_corr_init=zeros(nInit,nLag);
lagged_corr_detrend_init=zeros(nInit,nLag);

for target_day=1:nInit
for lag=0:nLag-1

corr_tmp=lagged_corr(target_day,lag+1);
corr_detrend_tmp=lagged_corr_detrend(target_day,lag+1);

init_day=target_day-lag;
if(init_day<=0)
init_day=init_day+365;
end

lagged_corr_init(init_day,lag+1)=corr_tmp;
lagged_corr_detrend_init(init_day,lag+1)=corr_detrend_tmp;

end
end


