function do_DTM_CHM

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% function do_DTM_CHM lataa hyytialan datat ja tekee niistä CHM         
%                  
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%      

eval(['load hyytiala2'])
file=first_pulse_hyytiala;

x_lower=min(file(:,2))-1;
y_lower=min(file(:,1))-1;
x_up=max(file(:,2))+1;
y_up=max(file(:,1))+1;
x_size=ceil(x_up-x_lower)*2;
y_size=ceil(y_up-y_lower)*2;
malli=zeros(x_size,y_size,2);
malli(:,:,1)=500;
l=length(file);
for i=1:l
    x_coord=file(i,2)-x_lower;
    y_coord=file(i,1)-y_lower;
    z_coord=file(i,3);
    x_coord=fix(x_coord*2+0.5);
    y_coord=fix(y_coord*2+0.5);
    if x_coord<x_size+1 & y_coord<y_size+1 & x_coord>0 & y_coord>0
      if malli(x_coord,y_coord,1)>z_coord
         malli(x_coord,y_coord,1)=z_coord;
      end
      if malli(x_coord,y_coord,2)<z_coord
         malli(x_coord,y_coord,2)=z_coord;
     end   
    end
end

file=last_pulse_hyytiala;

l=length(file);
for i=1:l
    x_coord=file(i,2)-x_lower;
    y_coord=file(i,1)-y_lower;
    z_coord=file(i,3);
    x_coord=fix(x_coord*2+0.5);
    y_coord=fix(y_coord*2+0.5);
    if x_coord<x_size+1 & y_coord<y_size+1 & x_coord>0 & y_coord>0
      if malli(x_coord,y_coord,1)>z_coord
         malli(x_coord,y_coord,1)=z_coord;
      end
      if malli(x_coord,y_coord,2)<z_coord
         malli(x_coord,y_coord,2)=z_coord;
     end   
    end
end

[col,row,d]=size(malli);
a=zeros(col,row);
osa=malli(:,:,1);
imagesc(osa)
pause
[x,y]=size(osa);
for i=1:x
   for j=1:y
      if i>10 & j>10 & i<x-10 & j<y-10
       		taul=osa(i-10:i+10,j-10:j+10);
            a(i,j)=min(min(taul));
      end
      b=0;
      if i>3 & j>3 & i<x-3 & j<y-3
         taul=osa(i-3:i+3,j-3:j+3);
         b=min(min(taul));
      end
      if abs(a(i,j)-b)<1
         a(i,j)=b;
      end 
      if abs(osa(i,j)-a(i,j))<1
         a(i,j)=osa(i,j);
      end  
   end
end
ero=osa-a;
[x,y]=find(abs(ero)<1);
z=osa(abs(ero)<1);
[x1,y1]=meshgrid(1:col,1:row);
clear osa
dtm=griddata(x,y,z,x1,y1,'linear')';
imagesc(dtm)
colorbar
pause
[x,y]=find(malli(:,:,2)>0);
z=malli(malli(:,:,2)>0);
dsm=griddata(x,y,z,x1,y1,'linear')';
imagesc(dsm)
colorbar
pause
chm=dsm-dtm;
imagesc(chm)
colorbar
pause

mal=malli(:,:,2)-dtm;
clear malli
malli=mal;

save hyytiala_lasketut2 chm malli x_lower y_lower