Skip to content
Snippets Groups Projects
Commit bf617fd7 authored by Gilles Boulet's avatar Gilles Boulet
Browse files

main_test.m simplified

parent a8485fed
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
Showing
with 37 additions and 152 deletions
main_test.m
+ 37
− 152
View file @ bf617fd7
......@@ -36,9 +36,9 @@ zf0=csvread('H_ble.csv');
LAI0=csvread('LAI_ble.csv');
LAIdate=LAI0(:,1); % dates de mesure de l'indice foliaire
gLAIval=LAI0(:,2); % valeur de l'indice foliaire ces dates
LAIval=LAI0(:,3); % valeur de l'indice foliaire total ces dates
LAIdate=LAI0(:,1); % Dates
gLAIval=LAI0(:,2); % Grean LAI at those dates
LAIval=LAI0(:,3); % Total LAI at those dates
glai=interp1(LAIdate,gLAIval,doy,'pchip','extrap');
lai=interp1(LAIdate,LAIval,doy,'pchip','extrap');
......@@ -94,7 +94,9 @@ betastress=0.01;
end
% bounding / stress
% bounding (can be adusted according to pot and stress baselines)
% bounding / stress
rns(les<lesstress)=rnsstress(les<lesstress);
g(les<lesstress)=gstress(les<lesstress);
......@@ -184,36 +186,23 @@ hpot2=hspot2+hvpot2;
a=doy-floor(doy);
time=find(a>0.485-(1/24) & a<0.515-(1/24));
time2=find(a>0.485-(1/24) & a<0.515-(1/24) & leobs >-100 & hobs>-100 & rnobs >-100 & gobs> - 100);
time3=time2;
RMSE_LE = sqrt(mean((leobs(time2) - le(time2)).^2)) % RMSE
RMSE_LE2 = sqrt(mean((leobs(time2) - le2(time2)).^2)) % RMSE 2
% Rott Mean Square Error on Latent Heat LE, Sensible HEat H, Net radiation
% Rn and Soil Heat Flux G
RMSE_H = sqrt(mean((hobs(time2) - h(time2)).^2)) % RMSE
RMSE_H2 = sqrt(mean((hobs(time2) - h2(time2)).^2)) % RMSE 2
RMSE_LE = sqrt(mean((leobs(time2) - le(time2)).^2)) % RMSE Layer
RMSE_LE2 = sqrt(mean((leobs(time2) - le2(time2)).^2)) % RMSE Patch
RMSE_RN = sqrt(mean((rnobs(time2) - rn(time2)).^2)) % RMSE
RMSE_RN2 = sqrt(mean((rnobs(time2) - rn2(time2)).^2)) % RMSE 2
RMSE_H = sqrt(mean((hobs(time2) - h(time2)).^2)) % RMSE Layer
RMSE_H2 = sqrt(mean((hobs(time2) - h2(time2)).^2)) % RMSE Patch
RMSE_G = sqrt(mean((gobs(time2) - g(time2)).^2)) % RMSE
RMSE_G2 = sqrt(mean((gobs(time2) - g2(time2)).^2)) % RMSE 2
RMSE_RN = sqrt(mean((rnobs(time2) - rn(time2)).^2)) % RMSE Layer
RMSE_RN2 = sqrt(mean((rnobs(time2) - rn2(time2)).^2)) % RMSE Patch
figure(1)
plot(doy(time),ts(time),'bo',doy(time),ts2(time),'bs',doy(time),tspot(time),'go',doy(time),tsstress(time),'ro',doy(time),tspot2(time),'gs',doy(time),tsstress2(time),'rs',doy(time),ta(time)-273.15,':k')
legend('simul srie','simul parallle','potentiel srie','stress srie','potentiel parallle','stress parallle','air')
title('temprature sol')
RMSE_G = sqrt(mean((gobs(time2) - g(time2)).^2)) % RMSE Layer
RMSE_G2 = sqrt(mean((gobs(time2) - g2(time2)).^2)) % RMSE Patch
figure(2)
plot(doy(time),tv(time),'bo',doy(time),tv2(time),'bs',doy(time),tvpot(time),'go',doy(time),tvstress(time),'ro',doy(time),tvpot2(time),'gs',doy(time),tvstress2(time),'rs',doy(time),ta(time)-273.15,':k')
legend('simul srie','simul parallle','potentiel srie','stress srie','potentiel parallle','stress parallle','air')
title('temprature vgtation')
figure(3)
plot(doy(time),trad(time),'bo',doy(time),trad2(time),'bs',doy(time),tradpot(time),'go',doy(time),tradstress(time),'ro',doy(time),tradpot2(time),'gs',doy(time),tradstress2(time),'rs',doy(time),ta(time)-273.15,':k',doy(time),tsobs(time),'ok')
legend('simul srie','simul parallle','potentiel srie','stress srie','potentiel parallle','stress parallle','air')
title('temprature radiative')
figure(4)
figure(1)
plot(leobs(time2),le(time2),'bo',leobs(time2),le2(time2),'rs')
hold on
lsline
......@@ -221,128 +210,51 @@ axis([-100 500 -100 500]);
legend('Series model','Parallel model')
xlabel('Observed Latent Heat Flux at midday [W/m2]')
ylabel('Retrieved Latent Heat Flux at midday [W/m2]')
% figure(12)
% %plot(doy(time),lev(time)./levpot(time),'ko',doy(time),lev2(time)./levpot(time),'bo',...
% % doy(time),les(time)./lespot(time),'k*',doy(time),les2(time)./lespot(time),'b*')
% % plot(doy(time),lev(time)./levpot(time),'go',doy(time),lev(time)./levstress(time),'gs',...
% % doy(time),les(time)./lespot(time),'ro',doy(time),les(time)./lesstress(time),'rs')
% plot(doy(time),Hv(time)./Hvpot(time),'go',doy(time),Hv(time)./Hvstress(time),'gs',...
% doy(time),Hs(time)./Hspot(time),'ro',doy(time),Hs(time)./Hsstress(time),'rs')
% axis([0 180 -1 2]);
figure(5)
plot(doy(time2),leobs(time2)./lepot(time2),'rs',doy(time2),le(time2)./lepot(time2),'bo',doy(time2),le2(time2)./lepot2(time2),'g*',doy(time2),lai(time2)/2,'g-')
%plot(doy(time2),leobs(time2).*2./(lepot(time2)+lepot2(time2)),'rs',doy(time2),lereal(time2)./lepot(time2),'ro',doy(time2),lereal2(time2)./lepot2(time2),'r*',doy(time2),le(time2)./lepot(time2),'bo',doy(time2),le2(time2)./lepot2(time2),'g*')
%axis([0 180 -0.2 1.2]);
legend('observed','series','parallel','LAI/2')
xlabel('Date [DOY]')
ylabel('Evaporative Efficiency LE/LEpot [-], LAI/2 [-]')
figure(6)
%plot(1-leobs(time3)./lepot(time3),1-le(time3)./lepot(time3),'bo',1-leobs(time3)./lepot2(time3),1-le2(time3)./lepot2(time3),'rs')
scatter(1-leobs(time3)./lepot(time3),1-le(time3)./lepot(time3),lepot(time3)/2,'bo')
hold on
scatter(1-leobs(time3)./lepot2(time3),1-le2(time3)./lepot2(time3),lepot2(time3)/2,'rs')
%lsline
axis([-0.2 1.2 -0.2 1.2]);
legend('Series model','Parallel model')
xlabel('Observed surface stress at midday (-)')
ylabel('Simulated surface stress at midday (-)')
figure(2)
scatter(1-leobs(time2)./lepot(time2),1-le(time2)./lepot(time2),'bo')
hold on
scatter(1-leobs(time2)./lepot2(time2),1-le2(time2)./lepot2(time2),'rs')
axis([-0.2 1.2 -0.2 1.2]);
legend('Series model','Parallel model')
xlabel('Observed surface stress at midday (-)')
ylabel('Simulated surface stress at midday (-)')
figure(7)
figure(3)
plot(doy(time2),leobs(time2),'rs',doy(time2),le(time2),'bo',doy(time2),le2(time2),'g*')
axis([0 180 0 500]);
legend('observ','srie','parallle')
xlabel('Date [DOY]')
ylabel('LE [W/m2]')
save('bounded.mat','doy','time2','le','le2')
weight=0.2*ones(1,5);
t=doy(time3);
es=les(time3);
es2=les2(time3);
esp=lespot(time3);
esp2=lespot2(time3);
t=doy(time2);
es=les(time2);
es2=les2(time2);
esp=lespot(time2);
esp2=lespot2(time2);
esf=filter(weight,1,es);
es2f=filter(weight,1,es2);
espf=filter(weight,1,esp);
esp2f=filter(weight,1,esp2);
figure(8)
%plot(doy(time2),les(time2),'bo',doy(time2),les2(time2),'g*',dateobs+0.5,pluie*10,'k-',doy(time2),lesobs(time2),'r:',doy(time2),leobs(time2),'rs')
figure(4)
plot(t,esf,'bo:',t,es2f,'rs:',t,espf,'g:',t,esp2f,'g-')
axis([0 180 0 250]);
legend('Series model','Parallel model','Potential serie','potential parallel')
xlabel('Date [DOY]')
ylabel('Soil Evaporation, 5 days running average [W/m2]')
figure(12)
figure(5)
plot(t,esf./espf,'bo:',t,es2f./esp2f,'rs:')
axis([0 180 -0.1 1.1]);
legend('Series model','Parallel model')
xlabel('Date [DOY]')
ylabel('Soil Evaporation Efficiency, 5 days running average [-]')
% figure(9)
% plot(doy(time2),lev2(time2),'bo',doy(time2),lev(time2),'b-',doy(time2),levpot(time2),'g:')
% legend('Old Srface Energy Balance Model TSEB','New Surface Energy Balance Model TSEB2','potentiel')
% xlabel('Date [DOY]')
% ylabel('LE [W/m2]')
figure(9)
plot(doy(time2),(h(time2)-hpot(time2))./(hstress(time2)-hpot(time2)),'b-',doy(time2),(hobs(time2)-hpot(time2))./(hstress(time2)-hpot(time2)),'bo',...
doy(time2),1-le(time2)./lepot(time2),'r-',doy(time2),1-leobs(time2)./lepot(time2),'rs')
%axis([0 180 -0.2 1.2]);
legend('Stress H simul','Stress H observ','Stress LE simul','Stress LE observ')
xlabel('Date [DOY]')
ylabel('Stress [-]')
figure(10)
plot(doy(time2),hv2(time2),'b-',doy(time2),hvpot2(time2),'r-',doy(time2),hvstress2(time2),'g-')
axis([0 180 -200 500]);
figure(11)
plot(doy(time2),hv(time2),'b-',doy(time2),hvpot(time2),'r-',doy(time2),hvstress(time2),'g-')
axis([0 180 -200 500]);
figure(13)
plot((hobs(time2)-hpot(time2))./(hstress(time2)-hpot(time2)),(h(time2)-hpot(time2))./(hstress(time2)-hpot(time2)),'bo',...
(hobs(time2)-hpot(time2))./(hstress(time2)-hpot(time2)),(h2(time2)-hpot(time2))./(hstress(time2)-hpot(time2)),'rs')
hold on
lsline
axis([-0.2 1.2 -0.2 1.2]);
legend('Series model','Parallel model')
xlabel('Observed surface stress at midday (-)')
ylabel('Simulated surface stress at midday (-)')
% figure(13)
% plot(doy(time2),leobs(time2),'rs',doy(time2),lev(time2),'bo',doy(time2),lev2(time2),'g*',dateobs+0.5,pluie*10,'k:')
% axis([0 180 -200 500]);
% legend('observ total','vgtation srie','vgtation parallle','pluie')
% xlabel('Date [DOY]')
% ylabel('LE [W/m2], pluie (*0.1 mm)')
% figure(14)
% plot(lesobs(time2),les(time2),'bo',lesobs(time2),les2(time2),'g*')
% hold on
% lsline
% axis([-100 400 -100 400]);
% title('vaporation du sol')
% legend('srie','parallle')
% xlabel('Flux observ midi (W/m2)')
% ylabel('Flux simul midi (W/m2)')
figure(15)
figure(6)
plot(hobs(time2),h(time2),'bo',hobs(time2),h2(time2),'g*')
hold on
lsline
......@@ -351,7 +263,7 @@ legend('series','parallel')
xlabel('Observed Sensible Heat Flux at Midday [W/m2]')
ylabel('Simulated Sensible Heat Flux at Midday [W/m2]')
figure(16)
figure(7)
plot(rnobs(time2),rn(time2),'bo',rnobs(time2),rn2(time2),'g*')
hold on
lsline
......@@ -359,7 +271,8 @@ axis([-100 700 -100 700]);
legend('series','parallel')
xlabel('Observed Net Radiation at Midday [W/m2]')
ylabel('Simulated Net Radition at Midday [W/m2]')
figure(18)
figure(8)
plot(gobs(time2),g(time2),'bo',gobs(time2),g2(time2),'g*')
hold on
lsline
......@@ -368,33 +281,5 @@ legend('series','parallel')
xlabel('Observed Soil Heat Flux at Midday [W/m2]')
ylabel('Observed Soil Heat Flux at Midday [W/m2]')
figure(21)
plot(doy(time3),1-leobs(time3)./lepot(time3),'rs',doy(time3),(tsobs(time3)-tradpot(time3))./(tradstress(time3)-tradpot(time3)),'b-',...
doy(time3),(tsobs(time3)-tradpot2(time3))./(tradstress2(time3)-tradpot2(time3)),'g-')
axis([0 180 -0.2 1.2]);
legend('Stress observ','Stress tso-tsp/ts0-tsp serie','Stress tso-tsp/ts0-tsp parallel')
xlabel('Date [DOY]')
ylabel('Stress [-]')
figure(22)
plot(1-leobs(time2)./lepot(time2),(tsobs(time2)-tradpot(time2))./(tradstress(time2)-tradpot(time2)),'bo',...
1-leobs(time2)./lepot(time2),1-le(time2)./lepot(time2),'go')
hold on
lsline
axis([-0.2 1.2 -0.2 1.2]);
legend('Stress tso-tsp/ts0-tsp serie','Stress 1-LE/LEpot')
xlabel('Stress observ [-]')
ylabel('Stress simul [-]')
figure(17)
plot(1-leobs(time2)./lepot(time2),(hobs(time2)-hpot(time2))./(hstress(time2)-hpot(time2)),'rs')
hold on
lsline
axis([-0.2 1.2 -0.2 1.2]);
legend('Stress LE','Stress H')
xlabel('Stress [-]')
ylabel('Stress [-]')
\ No newline at end of file
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment