function [yfit,res,h]=plot_weather_data_and_fit(x,y,f,p)
	% x points (date)
	% y points (temperature in Celsius degrees)
	% f - handle to fit function f(x,p)
        % p - fit parameters	
	% h -axes handles

	yf=@(x) f(x,p); % handle to fit function with parameters set
	yfit=arrayfun(yf,x);
	res=y-yfit;

	h(2)=subplot(2,1,2);
	plot(x,res,'.');  % plotting residuals
	hold on;
	plot(x,x*0, '-'); % this will plot a line along x axis
	legend('residuals','');
	hold off;
	xlabel('Date'); 
	ylabel('Residuals (^oC)');
	datetick('x','yyyy/mmm'); % x to date format


	h(1)=subplot(2,1,1); % notice I save hanle to axis !
	plot(x,y,'.'); % plotting data points
	hold on;
	plot(x,yfit,'r-'); % plotting fitted values
	hold off;
	legend('data','fit');
	ylabel('Temp (^oC)');
	datetick('x','yyyy/mmm'); % x to date format


	linkaxes(h, 'x' ); % link two subplot in 'x' direction

end
