%
% First, prompt user OR load in file with transect endpoints
%

cp=computer;


efgo=exist('formgrid.out');
if efgo~=0;
  wtd=input('A grid exist as "formgrid.out". What to do? 1=overwite,0=quit: ');
  if wtd==1;

     if cp(1:3)=='SGI';
        !rm formgrid.out;

     else;

        !del formgrid.out;

     end;   
  else;
    return;
  end;
end;

%iorb=menu('How do you want to specify end points?',...

%   'interactive','batch','graphic');
disp('How do you want to specify end points?');
iorb=input('Enter 1 for interactive, 2 for batch, or 3 for graphic: ');

first=input('enter first station# (usually 1): ');

% interactive option
if iorb==1|iorb==3;
    ns=input('number of transects? ');
    nsta=input('number of stations per transect? ');
end;


% batch option
if iorb==2;

    disp(' !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!')
    disp(' Program assumes you have a "formgrid.ctl"')
    disp(' in this directory which has five values per')
    disp(' record: start lat & lon, stop lat & lon, and')
    disp(' and desired spacing in nautical miles')
    disp(' and one record per transect.')
    disp('  ');
    disp(' PRESS RETURN TO CONTINUE NOW.')
    disp(' !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!')
    pause

    load formgrid.ctl
    [ns,mm]=size(formgrid);
%    if mm<6;
%      ns=input('number of transects? ');
%      nsta=input('number of stations per transect? ');
%    end;
    in=formgrid;
    clear formgrid;
  end;

%
j=1;
fid=fopen('formgrid.out','w');

% one transect at a time
for j=1:ns;

%note that if it is the first transect one must ask for positions
%but subsequent ones are calculated assuming same spacing in the
% perpendicular direction
  if j==1;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% interactive
    if iorb==1;
      slat=input('starting latitude in DDMM.T: ');
      slon=input('starting longitude in DDMM.T: ');
      elat=input('ending latitude in DDMM.T: ');
      elon=input('ending longitude in DDMM.T: ');
      space=input('Enter spacing in nautical miles: ');

    elseif iorb==3;

      disp(['click start point for transect #',num2str(j)]);

      [slon,slat]=ginput(1);

      disp(['click end point for transect #',num2str(j)]);

      [elon,elat]=ginput(1);

      [distkm,bear]=dist(slat,slon,elat,elon);

      % calculate spacing (miles) based on number of stations wanted

      space=distkm/nsta*60/111.209;
      disp(['Spacing is ',num2str(space),' miles']);
      space=input('Enter desired spacing in miles: ');
    elseif iorb==2;
      slat=in(j,1);
      slon=in(j,2);
      elat=in(j,3);
      elon=in(j,4);
      space=in(j,5);
    end;


% figure the angle perpendicular to this first section to be used
% on subsequent sections and also the lon/lat ratio for this region
    [distkm,bear]=dist(slat,slon,elat,elon);

% convert angle to radians counterclock wise from east
    bear=bear-90;
    bear=90-bear;
    bearrad=pi*(bear)/180;


% Convert these to decimal degrees
%
   [slat,slon]=ddmm2dec(slat,slon);
   [elat,elon]=ddmm2dec(elat,elon);

%   figure ratio after having converted to decdegree (below)
    [latd]=dist(slat,slon,slat-0.5,slon); 
    [lond]=dist(slat,slon,slat,slon+0.5); 
    llratio=latd/lond;

% figure lat step to next beginning section given 1st and bearing
    latstep=-1*space*sin(bearrad)/60;

% figure lon step by using the lon to lat distance ratio 
    lonstep=-1*(space*cos(bearrad)/60)*llratio*(-1);

    slat1=slat;
    slon1=slon;
    elat1=elat;
    elon1=elon;

% calculated from reference to first station
  elseif j>1;

    if iorb==1|iorb==3;
      slat=slat1+latstep*(j-1);
      slon=slon1-lonstep*(j-1); 
      elat=elat1+latstep*(j-1);
      elon=elon1-lonstep*(j-1); 
    end;
%   batch
    if iorb==2;
      slat=in(j,1);
      slon=in(j,2);
      elat=in(j,3);
      elon=in(j,4);
      space=in(j,5);
% Convert these to decimal degrees
%
   [slat,slon]=ddmm2dec(slat,slon);
   [elat,elon]=ddmm2dec(elat,elon);
    end;


  end;
    


% calculate how many stations ("milesint") should fit
  distkm=dist(slat,slon,elat,elon);
  milesint=floor(distkm/(1.81*space))+1;


% calculate the lat spacing and lon spacing
  latsp=(elat-slat)/(milesint-1);
  lonsp=(elon-slon)/(milesint-1);


% Increment lat & lon and write out stations
  ii=1;
  for ii=1:milesint;
    lat=slat+(ii-1)*latsp;
    lon=slon+(ii-1)*lonsp;
    z=w_depth(lon,lat);
% convert these to ddmm.t format
    [dmla,dmlo]=dd2dm(lat,lon);
    latdg=floor(dmla/100);latmm=dmla-latdg*100;
    londg=floor(dmlo/100);lonmm=dmlo-londg*100;

% write out a file to give to the bridge
    fprintf(fid,'%5.1f',first);
    fprintf(fid,'%8.0f %6.2f',latdg,latmm);
    fprintf(fid,'%8.0f %6.2f',londg,lonmm);
    fprintf(fid,'%10.3f %7.3f %3.0f',lat,lon,j);
    fprintf(fid,'%8.0f\n',z);
    first=first+1;
% go do another station
  end;


% go start another section
end;


fclose(fid);
