eha_rpt.m

function eha_rpt( p, t, s, ttle )
%EHA_RPT - Produces the final figure and results for the Eden-Hazel model
%                
% Syntax: eha_rpt( p, t, s, ttle )
%
%   p = parameters of the model
%   t = measured time
%   s = measured drawdown
%   ttle = Title of the figure 
%
% Description:
%   Produces the final figure and results for Eden-Hazel model.
%
% See also: eha_dmo, eha_pre, eha_gss, eha_dim
%


global EHA_MATRIX
if(  isempty(EHA_MATRIX) ) % NE MARCHE PAS - TROUVER UNE SOLUTION !
    disp(' ERROR: eha_rpt: You must run EHA_PRE before using the eha model')
    return;
end
d=EHA_MATRIX;

if(nargin==3) % Default value for d if not given by user
  ttle=' ';
end

% The p vector must be in line
if size(p,1) == 1        % It is a line vector then no problem 
elseif size(p,2) == 1    % It is a row vector then we have to transpose it 
  p=p';
else                     % It is a matrix, it cannot be treated
  error('p vector cannot be a matrix')
end

% Identification of the different pumping periods
np=d(end,2);  
for i=1:np
 j=find(d(:,2)==i);
 hl(i)=d(j(1),1);
 sl(i)=d(j(1),2);
 ql(i)=d(j(1),3);
end
np=d(end,2);

for i=1:np-1
 h(2*(i-1)+1)=hl(i);
 h(2*i)=hl(i+1);
 st(2*(i-1)+1)=sl(i);
 st(2*i)=sl(i);
 qi(2*(i-1)+1)=ql(i);
 qi(2*i)=ql(i);
end
h(2*(np-1)+1)=hl(np);
h(2*np)=d(end,1);
st(2*(np-1)+1)=sl(np);
st(2*np)=d(end,2);
qi(2*(np-1)+1)=ql(np);
qi(2*np)=d(end,3);

% Construction of the calculated Hn et drawdown
dplot=[h',st',qi'];  
sc=eha_cmp(p,dplot);

% Calculation of the hydraulic parameters
T=0.183/p(1);  % Transmissivity

% Calculation of the quadratic head losses coefficient C
A=p(2:end);      % Ai
X=qi(1:2:end)';  % X = flow rates = Qi
Y=(A'./X);         % Y = Ai/Qi
r=lin_fit(X,Y);   % Linear least-square Y=BX+A
B=r(1);
A=r(2);

% Plot the characteristic curve
clf
hax2=axes('position',[0.6,0.1,0.3,0.3]);
plot(X,Y,'o',X,B.*X+A)
xlabel('Qi')
ylabel('Ai/Qi')

% Select the area for the plot
xo=0.1; yo=0.5;
dy=0.4; dx=0.8;
hax=axes('position',[xo,yo,dx,dy]);
plot(d(:,1),s,'o',h,sc,'-')
xl=get(hax,'XLim');
yl=get(hax,'YLim');
legend('Drawdown','Model','Location','Northwest')
xlabel('Hn')
ylabel('Drawdown in m')
title(ttle)

% Coordinate transformation to position the legend
ax=(xl(2)-xl(1))./dx;
ay=(yl(2)-yl(1))./dy;

% Text of the legend
ht=text(ax*(0.1-xo)+xl(1),ay*(-0.1)+yl(1),'Eden-Hazel Model');
set(ht,'VerticalAlignment','top')
set(ht,'HorizontalAlignment','left')

legend1=sprintf('Fitting parameters:\n slope a: %0.2g m\n\nHydraulic parameters:\n Transmissivity T: %2.1e m^2/s\n Quadratic head losses B: %0.2g\n Constant A: %0.2g',p(1),T,B,A);
ht=text(ax*(0.1-xo)+xl(1),ay*(-0.18)+yl(1),legend1);
set(ht,'VerticalAlignment','top')
set(ht,'HorizontalAlignment','left')
 
copyright=hycoop;
hc=text(ax*(0.1-xo)+xl(1),ay*(-0.48)+yl(1),copyright);
set(hc,'FontSize',8)