model=1
%model=0;
load m2.matlab
%
if(model==1),
load model.m2
model(:,1)=model(:,1)*100;
model(:,2)=model(:,2)*100;
model(:,3)=90-model(:,3);
m2=m2*NaN;
m2(1,:)=model(1,:);
m2(3,:)=model(2,:);
m2(14,:)=model(3,:);
m2(21,:)=model(4,:);
m2(12,:)=model(5,:);
m2(16,:)=model(6,:);
m2(19,:)=model(7,:);
m2(24,:)=model(8,:);
m2(6,:)=model(9,:);
end

!rm rose.met
[m,n]=size(m2);
is=[1 3 6 12 14 16 19 21 24];   %surface 4-8 m
im=[2 4 7 13 15 17 20 22 25];   % mid 18-28 m
inb=[5 8 9  23];        %  33-60 m        % near bot 
ib=[10 11 18];           % > 60
ibb=[2  5  9 11 13 15  20  23]; % 10 m or less from  bottom
%
%  get coastline and 40 m isobath
%  set up plot
%lat_mid=42.3;
%lon_mid=-70.4;
%dlon=1.4;   
lat_mid=42.22;
lon_mid=-70.55;
dlon=1.1;
%lat_mid=42.22;
%lon_mid=-70.5;
%dlon=1.2;
lat_fac=cos(lat_mid*pi/180);
%
boxx=[-.5 .5 .5 -.5 -.5];
boxy=[-.5 -.5 .5 .5 -.5];
%
dlat=dlon*lat_fac;
% factor for roses
%
yfac=4*12.42*3600/(2*pi)*1.e-5/111.;   % [cm/s] * [s/day] * [km/cm] * [deg/km];
xfac=yfac/lat_fac;

clg
hold off
axis('square')
xmin=lon_mid-.5*dlon;
xmax=lon_mid+.5*dlon;
ymin=lat_mid-.5*dlat;
ymax=lat_mid+.5*dlat;
xx = [xmin xmax xmax xmin; xmax xmax xmin xmin];  % x coords of box
yy = [ymin ymin ymax ymax; ymin ymax ymax ymin];  % y coords of box

axis([xmin xmax ymin ymax])          % set axis limits
hold on
%
% str2='   4-8 m';
   str2=' 18-28 m';
% str2=' <= 10 mab';
for i=is,
% plot box
   xr=boxx*xfac+lonlat(i,1);
   yr=boxy*yfac+lonlat(i,2);
   plot(xr,yr,'r-')
% plot ellipse
   [z1]=tideell(m2(i,:));
   xr=real(z1)*xfac+lonlat(i,1);
   yr=imag(z1)*yfac+lonlat(i,2);
   plot(xr,yr,'g-')
%
end
if (exist('rose.met')~=2)
  meta rose
else
  meta
end
hold off
