function retstr = skyline_in_v3(ins)

% this function implements the topographic-shielding calculator.
%
% Syntax: out = skyline_in_v3(ins)
% 
% the input structure ins has the following fields, which all ought to be strings:
%
% ins.str_strike -- strike of surface in degrees
% ins.str_dip -- dip of surface in degrees
% ins.str_az -- space-separated string of azimuths
% ins.str_el -- space-separated string of elevations
%
% calls skyline.m to do the calculation
%
% Written by Greg Balco -- UW Cosmogenic Nuclide Lab
% balcs@bgc.org
% Originally written March, 2006
% Part of the CRONUS-Earth online calculators: 
%      http://hess.ess.washington.edu/math
% 
% Updated 2018
%
% Copyright 2001-2007, University of Washington
% subsequently, Greg Balco
% All rights reserved
% Developed in part with funding from the National Science Foundation.
%
% This program is free software; you can redistribute it and/or modify
% it under the terms of the GNU General Public License, version 2,
% as published by the Free Software Foundation (www.fsf.org).
%

% ------------------------------------------------

% What's the version?

ver = '2.0';

% Load constants

load /var/www/html/math/v3/consts_v3;

% Get in correct directory

if consts.isLocal == 0;
    % Case running on Linux web server
    cd /var/www/html/math/skyline
else
    % Case on GB's Mac, no cd but declare globals for diagnostics
    
end;

% ------------------------------------------------

% INPUT DATA CHECKING....

% Decide if anything was submitted
% These are booleans
isStrike =  ~(length(deblank(ins.str_strike)) == 0);
isDip =  ~(length(deblank(ins.str_dip)) == 0);
isAz = ~(length(deblank(ins.str_az)) == 0);
isEl = ~(length(deblank(ins.str_el)) == 0);

% check proper combinations

if ~((isStrike & isDip) == (isStrike | isDip));
    retstr = dump_error_HTML('skyline_in_v3.m: must enter both strike and dip');
    return;
end;


if ~((isAz & isEl) == (isAz | isEl));
    retstr = dump_error_HTML('skyline_in_v3.m: must enter both azimuths and elevations');
    return;
end;

% check for proper strike and dip values; do assignments

if isStrike;
    % Process whatever was entered for strike
    parsed_text = textscan(ins.str_strike,'%s');
    parsed_text = parsed_text{1}; % Peculiarity of textscan
    if length(parsed_text) > 1;
        retstr = dump_error_HTML('skyline_in_v3.m: strike value parses to multiple entries. Only one number please.');
        return;
    end;
    if isempty(regexp(parsed_text{1},'[^\d.eE+-]','once'));
		% pass, do assignment
		strike = str2double(parsed_text{1});
    else
		% fail
        retstr = dump_error_HTML('skyline_in_v3.m: unexpected characters in strike value -- use numbers only');
        return;
	end;
else
    strike = 0;
end;

clear parsed_text;
if isDip;
    % Process whatever was entered for dip
    parsed_text = textscan(ins.str_dip,'%s');
    parsed_text = parsed_text{1}; % Peculiarity of textscan
    if length(parsed_text) > 1;
        retstr = dump_error_HTML('skyline_in_v3.m: dip value parses to multiple entries. Only one number please.');
        return;
    end;
    if isempty(regexp(parsed_text{1},'[^\d.eE+-]','once'));
		% pass, do assignment
		dip = str2double(parsed_text{1});
    else
		% fail
        retstr = dump_error_HTML('skyline_in_v3.m: unexpected characters in dip value -- use numbers only');
        return;
	end;
else
    dip = 0;
end;

% Check for proper azimuth and elevation values; do assignments

az = []; el = []; % Allocate

if isAz;
    % Process whatever was entered for azimuth
    parsed_text = textscan(ins.str_az,'%s');
    parsed_text = parsed_text{1}; % Peculiarity of textscan
    for a = 1:length(parsed_text);
        if isempty(regexp(parsed_text{a},'[^\d.eE+-]','once'));
            % pass, do assignment
            az(end+1) = str2double(parsed_text{a});
        else
            % fail
            retstr = dump_error_HTML(['skyline_in_v3.m: unexpected characters in azimuth string, element ' int2str(a)]);
            return;
        end;
    end;
    % The result of that should be a vector of azimuths or an empty matrix.
    % Check range
    if any(az < 0) || any(az > 360);
        retstr = dump_error_HTML('skyline_in_v3.m: azimuths out of 0-360 range');
        return;
    end; 
end;

if isEl;
    % Process whatever was entered for azimuth
    parsed_text = textscan(ins.str_el,'%s');
    parsed_text = parsed_text{1}; % Peculiarity of textscan
    for a = 1:length(parsed_text);
        if isempty(regexp(parsed_text{a},'[^\d.eE+-]','once'));
            % pass, do assignment
            el(end+1) = str2double(parsed_text{a});
        else
            % fail
            retstr = dump_error_HTML(['skyline_in_v3.m: unexpected characters in elevation string, element ' int2str(a)]);
            return;
        end;
    end;
    % The result of that should be a vector of azimuths or an empty matrix.
    % Check range
    if any(el < 0) || any(el > 90);
        retstr = dump_error_HTML('skyline_in_v3.m: elevations out of 0-90 range');
        return;
    end; 
end;

if length(az) ~= length(el);
    retstr = dump_error_HTML('skyline_in_v3.m: azimuth and elevation arrays not the same size. Check for decimal points, commas, etc.');
    return;
end;

% Done checking input data. 

% Do shielding calculation

[sf,horiz,skyline_ver] = skyline(az,el,strike,dip);
% Do it for just surface and horizon separately for plot
[sf1,horiz1,skyline_ver] = skyline([],[],strike,dip);
[sf2,horiz2,skyline_ver] = skyline(az,el,0,0);

% --------- Use GMT to generate check plot ----------------------------

uniqueid = int2str(round(rand(1)*1e6));

fnamestub = ['skyline_' uniqueid];

if consts.isLocal == 0;
    % Assume running on hess.ess
    % everything goes on in the scratch directory.
    fname = ['/var/www/html/scratch/' fnamestub];
    gmtname = [fname '.gmt'];
    psname = [fname '.ps'];
    pngname = [fname '.png'];
else
    % Perhaps running on GB's laptop, refer to consts file
    fname = [consts.scratchdir fnamestub];
    gmtname = [fname '.gmt'];
    psname = [fname '.eps'];
    pngname = [fname '.png'];
end;
% open file for write
fid = fopen(gmtname,'w');

% Defs for GMT ver 5
pstr = 'gmt ';
wsep = ',';
formatstr = '--MAP_ANNOT_OFFSET_PRIMARY=0.1i --MAP_TICK_LENGTH_PRIMARY=-0.05i --MAP_FRAME_PEN=0.5p --FONT_ANNOT_PRIMARY=10p,2,black --FONT_LABEL=10p,3,black';

lxs = ['"Azimuth (degrees from N)"'];
lys = ['"Horizon angle (degrees)"'];


rstring = ['-R0/360/0/90'];
bstring = ['-Ba90g30:' lxs ':/a10g10:' lys ':WeSn']; % what is this?
jstring = '-JX5i/2i';

% Tell shell path to GMT if necessary, this is for GB's Mac
if consts.isLocal == 1;
    fprintf(fid,'%s\n\n','PATH=$PATH:/sw/bin');
    fprintf(fid,'%s\n\n','export PATH');
end;

% Make the GMT basemap
fprintf(fid,'%s\n','# ----------- Draw the plot axes ------------');

temp_string = ['gmt psbasemap ' bstring ' ' jstring ' ' rstring ' -K -P ' formatstr ' > ' psname];
fprintf(fid,'%s\n',temp_string);

% Make the surface horizon
temp_string = ['gmt psxy ' jstring ' ' rstring ' -W1p,255/0/0 -O -K -P << EOD >> ' psname];
fprintf(fid,'%s\n',temp_string);
fprintf(fid,'%0.5g %0.5g\n',[(1:360)' horiz1']');
fprintf(fid,'%s\n','EOD');

% Make the remote horizon
temp_string = ['gmt psxy ' jstring ' ' rstring ' -W1p,0/0/255 -O -K -P << EOD >> ' psname];
fprintf(fid,'%s\n',temp_string);
fprintf(fid,'%0.5g %0.5g\n',[(1:360)' horiz2']');
fprintf(fid,'%s\n','EOD');

% Make the complete horizon
temp_string = ['gmt psxy ' jstring ' ' rstring ' -W1p,0 -O -P << EOD >> ' psname];
fprintf(fid,'%s\n',temp_string);
fprintf(fid,'%0.5g %0.5g\n',[(1:360)' horiz']');
fprintf(fid,'%s\n','EOD');

% close the file
fclose(fid);

if consts.isLocal == 0;
    % Running on Linux server
    % Execute the script you just made
    system(['chmod a+x ' gmtname]);
    system(gmtname);

    % convert the ps file to a png using ImageMagick...
    system(['convert ' psname ' -trim -transparent white ' pngname]); 
    
else
    % On GB's Mac
    disp(['sh ' gmtname]);
    disp(['open ' psname]);
    % Running the convert command doesn't seem to work. 
end;

% return the name of the jpeg you just made
out.plotName = fnamestub;

% assemble other output data
out.sf = sf;
versions.wrapper = ver;
versions.skyline = skyline_ver;


% ------------------- WRITE TO LOG -----------------------------------

% Get IP address
this_ip = getenv('REMOTE_ADDR');

if consts.isLocal == 0;  
    % Write to the log file. 
    write_to_log_v3('skyline_in_v3.m',this_ip);
end;
    
% ------------------ DONE WRITING TO LOG -----------------------------

% -------------------- RETURN HTML ----------------------------

retstr = skylineToHTML(out,versions);

% -------------------------------- DONE ------------------------------
