% The path of residual file must be specified ..

% RESloc Function evaluation ....

function ZZ = NLLS_RESLOC_SIM(x,HHactual_1,HH_lag_1,EMPactual,sumHHzone,unemp,landv,landdRatio,landr,a_1,c,dimZone,dimHH,dimEMP,HH_lag);

%Initialization with the previous iteration values ....
eta=x(1,1); alpha=x(2,1); beta=x(3,1);
qq=x(4,1); rr=x(5,1); ss=x(6,1);
bnn=x(7:6+dimHH,1);

n=1;
%dimHH=1;

%Calculate Q 
Q=zeros(dimZone,n);
for j=1:dimZone;
        for k=1:dimEMP;
            Q(j,n)=Q(j,n)+(a_1(k,n).*EMPactual(j,k)); 
        end;
end;

%Calculate W_RESLOC
pipi=ones(dimZone,n); 
for i=1:dimZone;
        for nn=1:dimHH;
            pipi(i,n)=pipi(i,n).*((1+HH_lag(i,nn)./sumHHzone(i,1)).^bnn(nn,n)); 
        end;
  
        W_RESLOC(i,n)=landv(i,1).^qq(1,n)*(1+landdRatio(i,1)).^rr(1,n)*landr(i,1).^ss(1,n)*pipi(i,n);     
end;

%Calculate B
Binv=zeros(dimZone,n);
for j=1:dimZone  
        for i=1:dimZone
            Binv(j,n)=Binv(j,n)+(W_RESLOC(i,n).*c(i,j).^alpha(1,n).*exp(beta(1,n)*c(i,j)));
        end
        B(j,n)=1./Binv(j,n);
end

%Calculate modeled HH
HHsum=zeros(dimZone,n);
HH = zeros(dimZone,n);
for i=1:dimZone;
        for j=1:dimZone;
            HHsum(i,n)=HHsum(i,n)+Q(j,n).*B(j,n).*W_RESLOC(i,n).*c(i,j).^alpha(1,n).*exp(beta(1,n)*c(i,j));
        end;
        HH(i,n)=eta(1,n).*HHsum(i,n)+(1-eta(1,n)).*HH_lag_1(i,n);
        if(HH(i,n) < 0.001) HH(i,n) = 0.001;
        end;
end;


% Scaling the predicted values ..
sum_HHactual = sum(HHactual_1,1);
sum_HH = sum(HH,1);
HH = HH.*(sum_HHactual/sum_HH);


%Calculate residual
epsilon=zeros(dimZone,dimHH);
for i=1:dimZone;
        epsilon(i)= HHactual_1(i) - HH(i);
end;
csvwrite('C:\G-LUM Finale\Calibration\Output\epsilon_RESALLOC.dat',epsilon); % Outputpath must be specified

% Likelihood estimation
L = 0;
for i=1:dimZone;
    L = L + (HHactual_1(i).*log(HH(i,1)));
end;

% TELUM Objective .... 
ZZ = -L


