% look at average absolute errors of pixel values

x=avsread('../huber/phtm.fld'); xx=(1:128); [X Y]=meshgrid(xx,xx);

r1=and((x==12000),Y<60); r2=and((x==12000),Y>60); r3=(x==32000);
r4=or((x==22000),(x==26000));;
%r4=(X>70 & X<72 & Y<72 & Y>70);

xx=(-15:15); [X Y]=meshgrid(xx,xx); R2=X.^2+Y.^2;
sigma=2; Gf=1/sqrt(2*pi)/sigma*exp(-R2/2/sigma/sigma); Gf=Gf/sum(Gf(:));

y=avsread('../huber/phtmtbl.fld');
tt=y.*r1; t1=mean(tt(find(tt~=0)));
tt=y.*r2; t2=mean(tt(find(tt~=0)));
tt=y.*r3; t3=mean(tt(find(tt~=0)));
tt=y.*r4; t4=mean(tt(find(tt~=0)));

load ../myout5;
N=50; %z=avsread('../huber/hc/out11.fld');
z=myout5;
s1=zeros(1,N); s2=zeros(1,N); s3=zeros(1,N); s4=zeros(1,N); 
for i=1:50,
	%z(:,:,i)=filter2(Gf,myout4(:,:,i));
	ss=abs(z(:,:,i)-y).*r1; s1(i)=mean( ss(find(r1~=0)) );
        ss=abs(z(:,:,i)-y).*r2; s2(i)=mean( ss(find(r2~=0)) );
        ss=abs(z(:,:,i)-y).*r3; s3(i)=mean( ss(find(r3~=0)) );
        ss=abs(z(:,:,i)-y).*r4; s4(i)=mean( ss(find(r4~=0)) );
end
disp([mean(s1),(mean(s1))/t1,std(s1)*100]);
disp([mean(s2),(mean(s2))/t2,std(s2)*100]);
disp([mean(s3),(mean(s3))/t3,std(s3)*100]);
disp([mean(s4),(mean(s4))/t4,std(s4)*100]);

