+/****************************************************************************
+* SPLAT: A Terrain Analysis Program *
+* Copyright John A. Magliacane, KD2BD 1997-2004 *
+* Last update: 24-Jan-2004 *
+*****************************************************************************
+* *
+* This program is free software; you can redistribute it and/or modify it *
+* under the terms of the GNU General Public License as published by the *
+* Free Software Foundation; either version 2 of the License or any later *
+* version. *
+* *
+* This program is distributed in the hope that it will useful, but WITHOUT *
+* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or *
+* FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License *
+* for more details. *
+* *
+*****************************************************************************
+* *
+* Extensively modified by J. D. McDonald in Jan. 2004 to include *
+* the Longley-Rice propagation model using C++ code from NTIA/ITS. *
+* *
+* See: http://elbert.its.bldrdoc.gov/itm.html *
+* *
+*****************************************************************************
+ g++ -Wall -O3 -s -lm -lbz2 -fomit-frame-pointer itm.cpp splat.cpp -o splat
+*****************************************************************************/
+
+#include <stdio.h>
+#include <math.h>
+#include <stdlib.h>
+#include <string.h>
+#include <ctype.h>
+#include <bzlib.h>
+#include <unistd.h>
+#include "fontdata.h"
+#include "smallfont.h"
+
+#define GAMMA 2.5
+#define MAXSLOTS 9
+#define BZBUFFER 65536
+
+#if MAXSLOTS==4
+#define ARRAYSIZE 4950
+#endif
+
+#if MAXSLOTS==9
+#define ARRAYSIZE 10870
+#endif
+
+#if MAXSLOTS==16
+#define ARRAYSIZE 19240
+#endif
+
+#if MAXSLOTS==25
+#define ARRAYSIZE 30025
+#endif
+
+char string[255], sdf_path[255], opened=0, *splat_version={"1.1.0"};
+
+double TWOPI=6.283185307179586, HALFPI=1.570796326794896,
+ PI=3.141592653589793, deg2rad=1.74532925199e-02,
+ EARTHRADIUS=20902230.97, METERS_PER_MILE=1609.344,
+ METERS_PER_FOOT=0.3048, earthradius, max_range=0.0;
+
+int min_north=0, max_north=0, min_west=0, max_west=0,
+ max_elevation=0, min_elevation=0, bzerror;
+
+struct site { double lat;
+ double lon;
+ double alt;
+ char name[50];
+ };
+
+struct { float lat[ARRAYSIZE];
+ float lon[ARRAYSIZE];
+ float elevation[ARRAYSIZE];
+ float distance[ARRAYSIZE];
+ int length;
+ } path;
+
+struct { int min_north;
+ int max_north;
+ int min_west;
+ int max_west;
+ int max_el;
+ int min_el;
+ short data[1200][1200];
+ unsigned char mask[1200][1200];
+ } dem[MAXSLOTS];
+
+struct {
+ double eps_dielect;
+ double sgm_conductivity;
+ double eno_ns_surfref;
+ double frq_mhz;
+ double conf;
+ double rel;
+ int radio_climate;
+ int pol;
+ } LR;
+
+double elev_l[ARRAYSIZE+10];
+
+void point_to_point(double elev[], double tht_m, double rht_m,
+ double eps_dielect, double sgm_conductivity, double eno_ns_surfref,
+ double frq_mhz, int radio_climate, int pol, double conf,
+ double rel, double &dbloss, char *strmode, int &errnum);
+
+double arccos(double x, double y)
+{
+ /* This function implements the arc cosine function,
+ returning a value between 0 and TWOPI. */
+
+ double result=0.0;
+
+ if (y>0.0)
+ result=acos(x/y);
+
+ if (y<0.0)
+ result=PI+acos(x/y);
+
+ return result;
+}
+
+char *dec2dms(double decimal)
+{
+ /* Converts decimal degrees to degrees, minutes, seconds,
+ (DMS) and returns the result as a character string. */
+
+ int degrees, minutes, seconds;
+ double a, b, c, d;
+
+ a=floor(decimal);
+ b=60.0*(decimal-a);
+ c=floor(b);
+ d=60.0*(b-c);
+
+ degrees=(int)a;
+ minutes=(int)c;
+ seconds=(int)d;
+
+ if (seconds<0)
+ seconds=0;
+
+ if (seconds>59)
+ seconds=59;
+
+ string[0]=0;
+ sprintf(string,"%d%c %d\' %d\"", degrees, 176, minutes, seconds);
+ return (string);
+}
+
+int OrMask(double lat, double lon, int value)
+{
+ /* Lines, text, markings, and coverage areas are stored in a
+ mask that is combined with topology data when topographic
+ maps are generated by SPLAT!. This function sets bits in
+ the mask based on the latitude and longitude of the area
+ pointed to. */
+
+ int x, y, indx, minlat, minlon;
+ char found;
+
+ minlat=(int)floor(lat);
+ minlon=(int)floor(lon);
+
+ for (indx=0, found=0; indx<MAXSLOTS && found==0;)
+ if (minlat==dem[indx].min_north && minlon==dem[indx].min_west)
+ found=1;
+ else
+ indx++;
+
+ if (found)
+ {
+ x=(int)(1199.0*(lat-floor(lat)));
+ y=(int)(1199.0*(lon-floor(lon)));
+
+ dem[indx].mask[x][y]|=value;
+ return (dem[indx].mask[x][y]);
+ }
+ else
+ return -1;
+}
+
+int GetMask(double lat, double lon)
+{
+ /* This function returns the mask bits based on the latitude
+ and longitude given. */
+
+ return (OrMask(lat,lon,0));
+}
+
+double GetElevation(struct site location)
+{
+ /* This function returns the elevation (in feet) of any location
+ represented by the digital elevation model data in memory.
+ Function returns -5000.0 for locations not found in memory. */
+
+ char found;
+ int x, y, indx, minlat, minlon;
+ double elevation;
+
+ elevation=-5000.0;
+
+ minlat=(int)floor(location.lat);
+ minlon=(int)floor(location.lon);
+
+ x=(int)(1199.0*(location.lat-floor(location.lat)));
+ y=(int)(1199.0*(location.lon-floor(location.lon)));
+
+ for (indx=0, found=0; indx<MAXSLOTS && found==0; indx++)
+ {
+ if (minlat==dem[indx].min_north && minlon==dem[indx].min_west)
+ {
+ elevation=3.28084*dem[indx].data[x][y];
+ found=1;
+ }
+ }
+
+ return elevation;
+}
+
+double Distance(struct site site1, struct site site2)
+{
+ /* This function returns the great circle distance
+ in miles between any two site locations. */
+
+ double lat1, lon1, lat2, lon2, distance;
+
+ lat1=site1.lat*deg2rad;
+ lon1=site1.lon*deg2rad;
+ lat2=site2.lat*deg2rad;
+ lon2=site2.lon*deg2rad;
+
+ distance=3959.0*acos(sin(lat1)*sin(lat2)+cos(lat1)*cos(lat2)*cos((lon1)-(lon2)));
+
+ return distance;
+}
+
+double Azimuth(struct site source, struct site destination)
+{
+ /* This function returns the azimuth (in degrees) to the
+ destination as seen from the location of the source. */
+
+ double dest_lat, dest_lon, src_lat, src_lon,
+ beta, azimuth, diff, num, den, fraction;
+
+ dest_lat=destination.lat*deg2rad;
+ dest_lon=destination.lon*deg2rad;
+
+ src_lat=source.lat*deg2rad;
+ src_lon=source.lon*deg2rad;
+
+ /* Calculate Surface Distance */
+
+ beta=acos(sin(src_lat)*sin(dest_lat)+cos(src_lat)*cos(dest_lat)*cos(src_lon-dest_lon));
+
+ /* Calculate Azimuth */
+
+ num=sin(dest_lat)-(sin(src_lat)*cos(beta));
+ den=cos(src_lat)*sin(beta);
+ fraction=num/den;
+
+ /* Trap potential problems in acos() due to rounding */
+
+ if (fraction>=1.0)
+ fraction=1.0;
+
+ if (fraction<=-1.0)
+ fraction=-1.0;
+
+ /* Calculate azimuth */
+
+ azimuth=acos(fraction);
+
+ /* Reference it to True North */
+
+ diff=dest_lon-src_lon;
+
+ if (diff<=-PI)
+ diff+=TWOPI;
+
+ if (diff>=PI)
+ diff-=TWOPI;
+
+ if (diff>0.0)
+ azimuth=TWOPI-azimuth;
+
+ return (azimuth/deg2rad);
+}
+
+double ElevationAngle(struct site local, struct site remote)
+{
+ /* This function returns the angle of elevation (in degrees)
+ of the remote location as seen from the local site.
+ A positive result represents an angle of elevation (uptilt),
+ while a negative result represents an angle of depression
+ (downtilt), as referenced to a normal to the center of
+ the earth. */
+
+ register double a, b, dx;
+
+ a=GetElevation(remote)+remote.alt+earthradius;
+ b=GetElevation(local)+local.alt+earthradius;
+
+ dx=5280.0*Distance(local,remote);
+
+ /* Apply the Law of Cosines */
+
+ return ((180.0*(acos(((b*b)+(dx*dx)-(a*a))/(2.0*b*dx)))/PI)-90.0);
+}
+
+void ReadPath(struct site source, struct site destination)
+{
+ /* This function generates a sequence of latitude and
+ longitude positions between a source location and
+ a destination along a great circle path, and stores
+ elevation and distance information for points along
+ that path in the "path" structure for later use. */
+
+ double azimuth, distance, lat1, lon1, beta,
+ den, num, lat2, lon2, total_distance;
+ int x1, y1, c;
+ struct site tempsite;
+
+ c=0;
+
+ lat1=source.lat*deg2rad;
+ lon1=source.lon*deg2rad;
+ azimuth=Azimuth(source,destination)*deg2rad;
+ total_distance=Distance(source,destination);
+
+ for (distance=0; distance<=total_distance; distance+=0.04)
+ {
+ beta=distance/3959.0;
+ lat2=asin(sin(lat1)*cos(beta)+cos(azimuth)*sin(beta)*cos(lat1));
+ num=cos(beta)-(sin(lat1)*sin(lat2));
+ den=cos(lat1)*cos(lat2);
+
+ if (azimuth==0.0 && (beta>HALFPI-lat1))
+ lon2=lon1+PI;
+
+ else if (azimuth==HALFPI && (beta>HALFPI+lat1))
+ lon2=lon1+PI;
+
+ else if (fabs(num/den)>1.0)
+ lon2=lon1;
+
+ else
+ {
+ if ((PI-azimuth)>=0.0)
+ lon2=lon1-arccos(num,den);
+ else
+ lon2=lon1+arccos(num,den);
+ }
+
+ while (lon2<0.0)
+ lon2+=TWOPI;
+
+ while (lon2>TWOPI)
+ lon2-=TWOPI;
+
+ lat2=lat2/deg2rad;
+ lon2=lon2/deg2rad;
+
+ x1=(int)(1199.0*(lat2-floor(lat2)));
+ y1=(int)(1199.0*(lon2-floor(lon2)));
+
+ path.lat[c]=lat2;
+ path.lon[c]=lon2;
+ tempsite.lat=lat2;
+ tempsite.lon=lon2;
+ path.elevation[c]=GetElevation(tempsite);
+ path.distance[c]=distance;
+ c++;
+ }
+
+ /* Make sure exact destination point is recorded at path.length-1 */
+
+ x1=(int)(1199.0*(destination.lat-floor(destination.lat)));
+ y1=(int)(1199.0*(destination.lon-floor(destination.lon)));
+
+ path.lat[c]=destination.lat;
+ path.lon[c]=destination.lon;
+ path.elevation[c]=GetElevation(destination);
+ path.distance[c]=total_distance;
+ c++;
+
+ path.length=c;
+}
+
+double AverageTerrain(struct site source, double azimuthx, double start_distance, double end_distance)
+{
+ /* This function returns the average terrain calculated in
+ the direction of "azimuth" (degrees) between "start_distance"
+ and "end_distance" (miles) from the source location. If
+ the terrain is all water (non-critical error), -5000.0 is
+ returned. If not enough SDF data has been loaded into
+ memory to complete the survey (critical error), then
+ -9999.0 is returned. */
+
+ int c, samples, endpoint;
+ double beta, lat1, lon1, lat2, lon2, num, den, azimuth, terrain=0.0;
+ struct site destination;
+
+ lat1=source.lat*deg2rad;
+ lon1=source.lon*deg2rad;
+
+ /* Generate a path of elevations between the source
+ location and the remote location provided. */
+
+ beta=end_distance/3959.0;
+
+ azimuth=deg2rad*azimuthx;
+
+ lat2=asin(sin(lat1)*cos(beta)+cos(azimuth)*sin(beta)*cos(lat1));
+ num=cos(beta)-(sin(lat1)*sin(lat2));
+ den=cos(lat1)*cos(lat2);
+
+ if (azimuth==0.0 && (beta>HALFPI-lat1))
+ lon2=lon1+PI;
+
+ else if (azimuth==HALFPI && (beta>HALFPI+lat1))
+ lon2=lon1+PI;
+
+ else if (fabs(num/den)>1.0)
+ lon2=lon1;
+
+ else
+ {
+ if ((PI-azimuth)>=0.0)
+ lon2=lon1-arccos(num,den);
+ else
+ lon2=lon1+arccos(num,den);
+ }
+
+ while (lon2<0.0)
+ lon2+=TWOPI;
+
+ while (lon2>TWOPI)
+ lon2-=TWOPI;
+
+ lat2=lat2/deg2rad;
+ lon2=lon2/deg2rad;
+
+ destination.lat=lat2;
+ destination.lon=lon2;
+
+ /* If SDF data is missing for the endpoint of
+ the radial, then the average terrain cannot
+ be accurately calculated. Return -9999.0 */
+
+ if (GetElevation(destination)<-4999.0)
+ return (-9999.0);
+ else
+ {
+ ReadPath(source,destination);
+
+ endpoint=path.length;
+
+ /* Shrink the length of the radial if the
+ outermost portion is not over U.S. land. */
+
+ for (c=endpoint-1; c>=0 && path.elevation[c]==0.0; c--);
+
+ endpoint=c+1;
+
+ for (c=0, samples=0; c<endpoint; c++)
+ {
+ if (path.distance[c]>=start_distance)
+ {
+ terrain+=path.elevation[c];
+ samples++;
+ }
+ }
+
+ if (samples==0)
+ terrain=-5000.0; /* No land */
+ else
+ terrain=(terrain/(double)samples);
+
+ return terrain;
+ }
+}
+
+double haat(struct site antenna)
+{
+ /* This function returns the antenna's Height Above Average
+ Terrain (HAAT) based on FCC Part 73.313(d). If a critical
+ error occurs, such as a lack of SDF data to complete the
+ survey, -5000.0 is returned. */
+
+ int azi, c;
+ char error=0;
+ double terrain, avg_terrain, haat, sum=0.0;
+
+ /* Calculate the average terrain between 2 and 10 miles
+ from the antenna site at azimuths of 0, 45, 90, 135,
+ 180, 225, 270, and 315 degrees. */
+
+ for (c=0, azi=0; azi<=315 && error==0; azi+=45)
+ {
+ terrain=AverageTerrain(antenna, (double)azi, 2.0, 10.0);
+
+ if (terrain<-9998.0) /* SDF data is missing */
+ error=1;
+
+ if (terrain>-4999.0) /* It's land, not water */
+ {
+ sum+=terrain; /* Sum of averages */
+ c++;
+ }
+ }
+
+ if (error)
+ return -5000.0;
+ else
+ {
+ avg_terrain=(sum/(double)c);
+ haat=(antenna.alt+GetElevation(antenna))-avg_terrain;
+ return haat;
+ }
+}
+
+void PlaceMarker(struct site location)
+{
+ /* This function places text and marker data in the mask array
+ for illustration on topographic maps generated by SPLAT!.
+ By default, SPLAT! centers text information BELOW the marker,
+ but may move it above, to the left, or to the right of the
+ marker depending on how much room is available on the map,
+ or depending on whether the area is already occupied by
+ another marker or label. If no room or clear space is
+ available on the map to place the marker and its associated
+ text, then the marker and text are not written to the map. */
+
+ int a, b, c, byte;
+ char ok2print, occupied;
+ double x, y, lat, lon, textx=0.0, texty=0.0, xmin, xmax,
+ ymin, ymax, p1, p3, p6, p8, p12, p16, p24, label_length;
+
+ xmin=min_north;
+ xmax=max_north;
+ ymin=min_west;
+ ymax=max_west;
+ lat=location.lat;
+ lon=location.lon;
+
+ if (lat<xmax && lat>xmin && lon<ymax && lon>ymin)
+ {
+ p1=1.0/1200.0;
+ p3=3.0/1200.0;
+ p6=6.0/1200.0;
+ p8=8.0/1200.0;
+ p12=12.0/1200.0;
+ p16=16.0/1200.0;
+ p24=24.0/1200.0;
+ ok2print=0;
+ occupied=0;
+
+ /* Is Marker Position Clear Of Text Or Other Markers? */
+
+ for (x=lat-p3; (x<=xmax && x>=xmin && x<=lat+p3); x+=p1)
+ for (y=lon-p3; (y<=ymax && y>=ymin && y<=lon+p3); y+=p1)
+ occupied|=(GetMask(x,y)&2);
+
+ if (occupied==0)
+ {
+ /* Determine Where Text Can Be Positioned */
+
+ /* label_length=length in pixels.
+ Each character is 8 pixels wide. */
+
+ label_length=p1*(double)(strlen(location.name)<<3);
+
+ if (((lon+label_length)<=ymax) && (lon-label_length)>=ymin)
+ {
+ /* Default: Centered Text */
+
+ texty=lon+label_length/2.0;
+
+ if ((lat-p8)>=p16)
+ {
+ /* Position Text Below The Marker */
+
+ textx=lat-p8;
+
+ x=textx;
+ y=texty;
+
+ /* Is This Position Clear Of
+ Text Or Other Markers? */
+
+ for (a=0, occupied=0; a<16; a++)
+ {
+ for (b=0; b<(int)strlen(location.name); b++)
+ for (c=0; c<8; c++, y-=p1)
+ occupied|=(GetMask(x,y)&2);
+ x-=p1;
+ y=texty;
+ }
+
+ x=textx;
+ y=texty;
+
+ if (occupied==0)
+ ok2print=1;
+ }
+
+ else
+ {
+ /* Position Text Above The Marker */
+
+ textx=lat+p24;
+
+ x=textx;
+ y=texty;
+
+ /* Is This Position Clear Of
+ Text Or Other Markers? */
+
+ for (a=0, occupied=0; a<16; a++)
+ {
+ for (b=0; b<(int)strlen(location.name); b++)
+ for (c=0; c<8; c++, y-=p1)
+ occupied|=(GetMask(x,y)&2);
+ x-=p1;
+ y=texty;
+ }
+
+ x=textx;
+ y=texty;
+
+ if (occupied==0)
+ ok2print=1;
+ }
+ }
+
+ if (ok2print==0)
+ {
+ if ((lon-label_length)>=ymin)
+ {
+ /* Position Text To The
+ Right Of The Marker */
+
+ textx=lat+p6;
+ texty=lon-p12;
+
+ x=textx;
+ y=texty;
+
+ /* Is This Position Clear Of
+ Text Or Other Markers? */
+
+ for (a=0, occupied=0; a<16; a++)
+ {
+ for (b=0; b<(int)strlen(location.name); b++)
+ for (c=0; c<8; c++, y-=p1)
+ occupied|=(GetMask(x,y)&2);
+ x-=p1;
+ y=texty;
+ }
+
+ x=textx;
+ y=texty;
+
+ if (occupied==0)
+ ok2print=1;
+ }
+
+ else
+ {
+ /* Position Text To The
+ Left Of The Marker */
+
+ textx=lat+p6;
+ texty=lon+p8+(label_length);
+
+ x=textx;
+ y=texty;
+
+ /* Is This Position Clear Of
+ Text Or Other Markers? */
+
+ for (a=0, occupied=0; a<16; a++)
+ {
+ for (b=0; b<(int)strlen(location.name); b++)
+ for (c=0; c<8; c++, y-=p1)
+ occupied|=(GetMask(x,y)&2);
+ x-=p1;
+ y=texty;
+ }
+
+ x=textx;
+ y=texty;
+
+ if (occupied==0)
+ ok2print=1;
+ }
+ }
+
+ /* textx and texty contain the latitude and longitude
+ coordinates that describe the placement of the text
+ on the map. */
+
+ if (ok2print && textx!=0.0 && texty!=0.0)
+ {
+ /* Draw Text */
+
+ x=textx;
+ y=texty;
+
+ for (a=0; a<16 && ok2print; a++)
+ {
+ for (b=0; b<(int)strlen(location.name); b++)
+ {
+ byte=fontdata[16*(location.name[b])+a];
+
+ for (c=128; c>0; c=c>>1, y-=p1)
+ if (byte&c)
+ OrMask(x,y,2);
+ }
+ x-=p1;
+ y=texty;
+ }
+
+ /* Draw Square Marker Centered
+ On Location Specified */
+
+ for (x=lat-p3; (x<=xmax && x>=xmin && x<=lat+p3); x+=p1)
+ for (y=lon-p3; (y<=ymax && y>=ymin && y<=lon+p3); y+=p1)
+ OrMask(x,y,2);
+ }
+ }
+ }
+}
+
+double ReadBearing(char *input)
+{
+ /* This function takes numeric input in the form of a character
+ string, and returns an equivalent bearing in degrees as a
+ decimal number (double). The input may either be expressed
+ in decimal format (40.139722) or degree, minute, second
+ format (40 08 23). This function also safely handles
+ extra spaces found either leading, trailing, or
+ embedded within the numbers expressed in the
+ input string. Decimal seconds are permitted. */
+
+ double bearing=0.0;
+ char string[20];
+ int a, b, length, degrees, minutes;
+ double seconds;
+
+ /* Copy "input" to "string", and ignore any extra
+ spaces that might be present in the process. */
+
+ string[0]=0;
+ length=strlen(input);
+
+ for (a=0, b=0; a<length && a<18; a++)
+ {
+ if ((input[a]!=32 && input[a]!='\n') || (input[a]==32 && input[a+1]!=32 && b!=0))
+ {
+ string[b]=input[a];
+ b++;
+ }
+ }
+
+ string[b]=0;
+
+ /* Count number of spaces in the clean string. */
+
+ length=strlen(string);
+
+ for (a=0, b=0; a<length; a++)
+ if (string[a]==32)
+ b++;
+
+ if (b==0) /* Decimal Format (40.139722) */
+ sscanf(string,"%lf",&bearing);
+
+ if (b==2) /* Degree, Minute, Second Format (40 08 23) */
+ {
+ sscanf(string,"%d %d %lf",°rees, &minutes, &seconds);
+ bearing=(double)degrees+((double)minutes/60)+(seconds/3600);
+ }
+
+ /* Anything else returns a 0.0 */
+
+ if (bearing>360.0 || bearing<0.0)
+ bearing=0.0;
+
+ return bearing;
+}
+
+struct site LoadQTH(char *filename)
+{
+ /* This function reads SPLAT! .qth (site location) files.
+ The latitude and longitude may be expressed either in
+ decimal degrees, or in degree, minute, second format.
+ Antenna height is assumed to be expressed in feet above
+ ground level (AGL), unless followed by the letter 'M',
+ or 'm', or by the word "meters" or "Meters", in which
+ case meters is assumed, and is handled accordingly. */
+
+ int x;
+ char string[50], qthfile[255];
+ struct site tempsite;
+ FILE *fd=NULL;
+
+ for (x=0; filename[x]!='.' && filename[x]!=0 && x<250; x++)
+ qthfile[x]=filename[x];
+
+ qthfile[x]='.';
+ qthfile[x+1]='q';
+ qthfile[x+2]='t';
+ qthfile[x+3]='h';
+ qthfile[x+4]=0;
+
+ tempsite.lat=0.0;
+ tempsite.lon=0.0;
+ tempsite.alt=0.0;
+ tempsite.name[0]=0;
+
+ fd=fopen(qthfile,"r");
+
+ if (fd!=NULL)
+ {
+ /* Site Name */
+ fgets(string,49,fd);
+
+ /* Strip <CR> and/or <LF> from end of site name */
+
+ for (x=0; string[x]!=13 && string[x]!=10 && string[x]!=0; tempsite.name[x]=string[x], x++);
+
+ tempsite.name[x]=0;
+
+ /* Site Latitude */
+ fgets(string,49,fd);
+ tempsite.lat=ReadBearing(string);
+
+ /* Site Longitude */
+ fgets(string,49,fd);
+ tempsite.lon=ReadBearing(string);
+
+ /* Antenna Height */
+ fgets(string,49,fd);
+ fclose(fd);
+
+ /* Remove <CR> and/or <LF> from antenna height string */
+ for (x=0; string[x]!=13 && string[x]!=10 && string[x]!=0; x++);
+
+ string[x]=0;
+
+ /* Antenna height may either be in feet or meters.
+ If the letter 'M' or 'm' is discovered in
+ the string, then this is an indication that
+ the value given is expressed in meters, and
+ must be converted to feet before exiting. */
+
+ for (x=0; string[x]!='M' && string[x]!='m' && string[x]!=0 && x<48; x++);
+ if (string[x]=='M' || string[x]=='m')
+ {
+ string[x]=0;
+ sscanf(string,"%lf",&tempsite.alt);
+ tempsite.alt*=3.28084;
+ }
+
+ else
+ {
+ string[x]=0;
+ sscanf(string,"%lf",&tempsite.alt);
+ }
+ }
+
+ return tempsite;
+}
+
+int LoadSDF_SDF(char *name)
+{
+ /* This function reads uncompressed SPLAT Data Files (.sdf)
+ containing digital elevation model data into memory.
+ Elevation data, maximum and minimum elevations, and
+ quadrangle limits are stored in the first available
+ dem[] structure. */
+
+ int x, y, data, indx, minlat, minlon, maxlat, maxlon;
+ char found, free_slot=0, line[20], sdf_file[255], path_plus_name[255];
+ FILE *fd;
+
+ for (x=0; name[x]!='.' && name[x]!=0 && x<250; x++)
+ sdf_file[x]=name[x];
+
+ sdf_file[x]=0;
+
+ /* Parse filename for minimum latitude and longitude values */
+
+ sscanf(sdf_file,"%d:%d:%d:%d",&minlat,&maxlat,&minlon,&maxlon);
+
+ sdf_file[x]='.';
+ sdf_file[x+1]='s';
+ sdf_file[x+2]='d';
+ sdf_file[x+3]='f';
+ sdf_file[x+4]=0;
+
+ /* Is it already in memory? */
+
+ for (indx=0, found=0; indx<MAXSLOTS && found==0; indx++)
+ {
+ if (minlat!=0 && minlon!=0)
+ {
+ if (minlat==dem[indx].min_north && minlon==dem[indx].min_west && maxlat==dem[indx].max_north && maxlon==dem[indx].max_west)
+ found=1;
+ }
+ }
+
+ /* Is room available to load it? */
+
+ if (found==0)
+ {
+ for (indx=0, free_slot=0; indx<MAXSLOTS && free_slot==0; indx++)
+ if (dem[indx].max_north==0 && dem[indx].max_west==0)
+ free_slot=1;
+ }
+
+ indx--;
+
+ if (free_slot && found==0 && indx>=0 && indx<MAXSLOTS && minlat!=0 && minlon!=0)
+ {
+ strncpy(path_plus_name,sdf_path,255);
+ strncat(path_plus_name,sdf_file,255);
+
+ fd=fopen(path_plus_name,"rb");
+
+ if (fd!=NULL)
+ {
+ fprintf(stdout,"Loading \"%s\" into slot %d...",path_plus_name,indx+1);
+ fflush(stdout);
+
+ fgets(line,19,fd);
+ sscanf(line,"%d",&dem[indx].max_west);
+
+ fgets(line,19,fd);
+ sscanf(line,"%d",&dem[indx].min_north);
+
+ fgets(line,19,fd);
+ sscanf(line,"%d",&dem[indx].min_west);
+
+ fgets(line,19,fd);
+ sscanf(line,"%d",&dem[indx].max_north);
+
+ for (x=0; x<1200; x++)
+ for (y=0; y<1200; y++)
+ {
+ fgets(line,19,fd);
+ sscanf(line,"%d",&data);
+ dem[indx].data[x][y]=data;
+
+ if (data>dem[indx].max_el)
+ dem[indx].max_el=data;
+
+ if (dem[indx].min_el==0)
+ dem[indx].min_el=data;
+ else
+ {
+ if (data<dem[indx].min_el)
+ dem[indx].min_el=data;
+ }
+ }
+
+ fclose(fd);
+
+ if (min_elevation==0)
+ min_elevation=dem[indx].min_el;
+
+ else
+ {
+ if (dem[indx].min_el<min_elevation)
+ min_elevation=dem[indx].min_el;
+ }
+
+ if (dem[indx].max_el>max_elevation)
+ max_elevation=dem[indx].max_el;
+
+ if (dem[indx].max_north>max_north)
+ max_north=dem[indx].max_north;
+
+ if (dem[indx].max_west>max_west)
+ max_west=dem[indx].max_west;
+
+ if (min_north==0)
+ min_north=dem[indx].min_north;
+ else
+ {
+ if (dem[indx].min_north<min_north)
+ min_north=dem[indx].min_north;
+ }
+
+ if (min_west==0)
+ min_west=dem[indx].min_west;
+ else
+ {
+ if (dem[indx].min_west<min_west)
+ min_west=dem[indx].min_west;
+ }
+
+ fprintf(stdout," Done!\n");
+ fflush(stdout);
+ return 1;
+ }
+
+ else
+ return -1;
+ }
+
+ else
+ return 0;
+}
+
+char *BZfgets(BZFILE *bzfd, unsigned length)
+{
+ /* This function returns at most one less than 'length' number
+ of characters from a bz2 compressed file whose file descriptor
+ is pointed to by *bzfd. In operation, a buffer is filled with
+ uncompressed data (size = BZBUFFER), which is then parsed
+ and doled out as NULL terminated character strings every time
+ this function is invoked. A NULL string indicates an EOF
+ or error condition. */
+
+ static int x, y, nBuf;
+ static char buffer[BZBUFFER+1], output[BZBUFFER+1];
+ char done=0;
+
+ if (opened!=1 && bzerror==BZ_OK)
+ {
+ /* First time through. Initialize everything! */
+
+ x=0;
+ y=0;
+ nBuf=0;
+ opened=1;
+ output[0]=0;
+ }
+
+ do
+ {
+ if (x==nBuf && bzerror!=BZ_STREAM_END && bzerror==BZ_OK && opened)
+ {
+ /* Uncompress data into a static buffer */
+
+ nBuf=BZ2_bzRead(&bzerror, bzfd, buffer, BZBUFFER);
+ buffer[nBuf]=0;
+ x=0;
+ }
+
+ /* Build a string from buffer contents */
+
+ output[y]=buffer[x];
+
+ if (output[y]=='\n' || output[y]==0 || y==(int)length-1)
+ {
+ output[y+1]=0;
+ done=1;
+ y=0;
+ }
+
+ else
+ y++;
+ x++;
+
+ } while (done==0);
+
+ if (output[0]==0)
+ opened=0;
+
+ return (output);
+}
+
+int LoadSDF_BZ(char *name)
+{
+ /* This function reads .bz2 compressed SPLAT Data Files containing
+ digital elevation model data into memory. Elevation data,
+ maximum and minimum elevations, and quadrangle limits are
+ stored in the first available dem[] structure. */
+
+ int x, y, data, indx, minlat, minlon, maxlat, maxlon;
+ char found, free_slot=0, sdf_file[255], path_plus_name[255];
+ FILE *fd;
+ BZFILE *bzfd;
+
+ for (x=0; name[x]!='.' && name[x]!=0 && x<247; x++)
+ sdf_file[x]=name[x];
+
+ sdf_file[x]=0;
+
+ /* Parse sdf_file name for minimum latitude and longitude values */
+
+ sscanf(sdf_file,"%d:%d:%d:%d",&minlat,&maxlat,&minlon,&maxlon);
+
+ sdf_file[x]='.';
+ sdf_file[x+1]='s';
+ sdf_file[x+2]='d';
+ sdf_file[x+3]='f';
+ sdf_file[x+4]='.';
+ sdf_file[x+5]='b';
+ sdf_file[x+6]='z';
+ sdf_file[x+7]='2';
+ sdf_file[x+8]=0;
+
+ /* Is it already in memory? */
+
+ for (indx=0, found=0; indx<MAXSLOTS && found==0; indx++)
+ {
+ if (minlat!=0 && minlon!=0)
+ {
+ if (minlat==dem[indx].min_north && minlon==dem[indx].min_west && maxlat==dem[indx].max_north && maxlon==dem[indx].max_west)
+ found=1;
+ }
+ }
+
+ /* Is room available to load it? */
+
+ if (found==0)
+ {
+ for (indx=0, free_slot=0; indx<MAXSLOTS && free_slot==0; indx++)
+ if (dem[indx].max_north==0 && dem[indx].max_west==0)
+ free_slot=1;
+ }
+
+ indx--;
+
+ if (free_slot && found==0 && indx>=0 && indx<MAXSLOTS && minlat!=0 && minlon!=0)
+ {
+ strncpy(path_plus_name,sdf_path,255);
+ strncat(path_plus_name,sdf_file,255);
+
+ fd=fopen(path_plus_name,"rb");
+ bzfd=BZ2_bzReadOpen(&bzerror,fd,0,0,NULL,0);
+
+ if (fd!=NULL && bzerror==BZ_OK)
+ {
+ fprintf(stdout,"Loading \"%s\" into slot %d...",path_plus_name,indx+1);
+ fflush(stdout);
+
+ sscanf(BZfgets(bzfd,255),"%d",&dem[indx].max_west);
+ sscanf(BZfgets(bzfd,255),"%d",&dem[indx].min_north);
+ sscanf(BZfgets(bzfd,255),"%d",&dem[indx].min_west);
+ sscanf(BZfgets(bzfd,255),"%d",&dem[indx].max_north);
+
+ for (x=0; x<1200; x++)
+ for (y=0; y<1200; y++)
+ {
+ sscanf(BZfgets(bzfd,20),"%d",&data);
+ dem[indx].data[x][y]=data;
+
+ if (data>dem[indx].max_el)
+ dem[indx].max_el=data;
+
+ if (dem[indx].min_el==0)
+ dem[indx].min_el=data;
+ else
+ {
+ if (data<dem[indx].min_el)
+ dem[indx].min_el=data;
+ }
+ }
+
+ fclose(fd);
+
+ BZ2_bzReadClose(&bzerror,bzfd);
+
+ if (min_elevation==0)
+ min_elevation=dem[indx].min_el;
+
+ else
+ {
+ if (dem[indx].min_el<min_elevation)
+ min_elevation=dem[indx].min_el;
+ }
+
+ if (dem[indx].max_el>max_elevation)
+ max_elevation=dem[indx].max_el;
+
+ if (dem[indx].max_north>max_north)
+ max_north=dem[indx].max_north;
+
+ if (dem[indx].max_west>max_west)
+ max_west=dem[indx].max_west;
+
+ if (min_north==0)
+ min_north=dem[indx].min_north;
+ else
+ {
+ if (dem[indx].min_north<min_north)
+ min_north=dem[indx].min_north;
+ }
+
+ if (min_west==0)
+ min_west=dem[indx].min_west;
+ else
+ {
+ if (dem[indx].min_west<min_west)
+ min_west=dem[indx].min_west;
+ }
+
+ fprintf(stdout," Done!\n");
+ fflush(stdout);
+ return 1;
+ }
+ else
+ return -1;
+ }
+ else
+ return 0;
+}
+
+char LoadSDF(char *name)
+{
+ /* This function loads the requested SDF file from the filesystem.
+ It first tries to invoke the LoadSDF_SDF() function to load an
+ uncompressed SDF file (since uncompressed files load slightly
+ faster). If that attempt fails, then it tries to load a
+ compressed SDF file by invoking the LoadSDF_BZ() function.
+ If that fails, then we can assume that no elevation data
+ exists for the region requested, and that the region
+ requested must be entirely over water. */
+
+ int x, y, indx, minlat, minlon, maxlat, maxlon;
+ char found, free_slot=0;
+ int return_value=-1;
+
+ /* Try to load an uncompressed SDF first. */
+
+ return_value=LoadSDF_SDF(name);
+
+ /* If that fails, try loading a compressed SDF. */
+
+ if (return_value==0 || return_value==-1)
+ return_value=LoadSDF_BZ(name);
+
+ /* If neither format can be found, then assume the area is water. */
+
+ if (return_value==0 || return_value==-1)
+ {
+ /* Parse SDF name for minimum latitude and longitude values */
+
+ sscanf(name,"%d:%d:%d:%d",&minlat,&maxlat,&minlon,&maxlon);
+
+ /* Is it already in memory? */
+
+ for (indx=0, found=0; indx<MAXSLOTS && found==0; indx++)
+ {
+ if (minlat!=0 && minlon!=0)
+ {
+ if (minlat==dem[indx].min_north && minlon==dem[indx].min_west && maxlat==dem[indx].max_north && maxlon==dem[indx].max_west)
+ found=1;
+ }
+ }
+
+ /* Is room available to load it? */
+
+ if (found==0)
+ {
+ for (indx=0, free_slot=0; indx<MAXSLOTS && free_slot==0; indx++)
+ if (dem[indx].max_north==0 && dem[indx].max_west==0)
+ free_slot=1;
+ }
+
+ indx--;
+
+ if (free_slot && found==0 && indx>=0 && indx<MAXSLOTS && minlat!=0 && minlon!=0)
+ {
+ fprintf(stdout,"Region \"%s\" assumed as sea-level into slot %d...",name,indx+1);
+ fflush(stdout);
+
+ dem[indx].max_west=maxlon;
+ dem[indx].min_north=minlat;
+ dem[indx].min_west=minlon;
+ dem[indx].max_north=maxlat;
+
+ for (x=0; x<1200; x++)
+ for (y=0; y<1200; y++)
+ {
+ dem[indx].data[x][y]=0;
+
+ if (dem[indx].min_el>0)
+ dem[indx].min_el=0;
+ }
+
+ if (dem[indx].min_el<min_elevation)
+ min_elevation=dem[indx].min_el;
+
+ if (dem[indx].max_el>max_elevation)
+ max_elevation=dem[indx].max_el;
+
+ if (dem[indx].max_north>max_north)
+ max_north=dem[indx].max_north;
+
+ if (dem[indx].max_west>max_west)
+ max_west=dem[indx].max_west;
+
+ if (min_north==0)
+ min_north=dem[indx].min_north;
+ else
+ {
+ if (dem[indx].min_north<min_north)
+ min_north=dem[indx].min_north;
+ }
+
+ if (min_west==0)
+ min_west=dem[indx].min_west;
+ else
+ {
+ if (dem[indx].min_west<min_west)
+ min_west=dem[indx].min_west;
+ }
+
+ fprintf(stdout," Done!\n");
+ fflush(stdout);
+
+ return_value=1;
+ }
+ }
+
+ return return_value;
+}
+
+void LoadCities(char *filename)
+{
+ /* This function reads SPLAT! city/site files, and plots
+ the locations and names of the cities and site locations
+ read on topographic maps generated by SPLAT! */
+
+ int x, y, z;
+ char input[80], str[3][80];
+ struct site city_site;
+ FILE *fd=NULL;
+
+ fd=fopen(filename,"r");
+
+ if (fd!=NULL)
+ {
+ fgets(input,78,fd);
+
+ fprintf(stdout,"Reading \"%s\"... ",filename);
+ fflush(stdout);
+
+ while (fd!=NULL && feof(fd)==0)
+ {
+ /* Parse line for name, latitude, and longitude */
+
+ for (x=0, y=0, z=0; x<78 && input[x]!=0 && z<3; x++)
+ {
+ if (input[x]!=',' && y<78)
+ {
+ str[z][y]=input[x];
+ y++;
+ }
+
+ else
+ {
+ str[z][y]=0;
+ z++;
+ y=0;
+ }
+ }
+
+ strncpy(city_site.name,str[0],49);
+ city_site.lat=ReadBearing(str[1]);
+ city_site.lon=ReadBearing(str[2]);
+ city_site.alt=0.0;
+
+ PlaceMarker(city_site);
+
+ fgets(input,78,fd);
+ }
+
+ fclose(fd);
+ fprintf(stdout,"Done!\n");
+ fflush(stdout);
+ }
+ else
+ fprintf(stderr,"*** ERROR: \"%s\": not found!\n",filename);
+}
+
+void LoadBoundaries(char *filename)
+{
+ /* This function reads Cartographic Boundary Files available from
+ the U.S. Census Bureau, and plots the data contained in those
+ files on the PPM Map generated by SPLAT!. Such files contain
+ the coordinates that describe the boundaries of cities,
+ counties, and states. */
+
+ int x;
+ double lat0, lon0, lat1, lon1;
+ char string[80];
+ struct site source, destination;
+ FILE *fd=NULL;
+
+ fd=fopen(filename,"r");
+
+ if (fd!=NULL)
+ {
+ fgets(string,78,fd);
+
+ fprintf(stdout,"Reading \"%s\"... ",filename);
+ fflush(stdout);
+
+ do
+ {
+ fgets(string,78,fd);
+ sscanf(string,"%lf %lf", &lon0, &lat0);
+ fgets(string,78,fd);
+
+ do
+ {
+ sscanf(string,"%lf %lf", &lon1, &lat1);
+
+ lon0=fabs(lon0);
+ lon1=fabs(lon1);
+
+ source.lat=lat0;
+ source.lon=lon0;
+ destination.lat=lat1;
+ destination.lon=lon1;
+
+ ReadPath(source,destination);
+
+ for (x=0; x<path.length; x++)
+ OrMask(path.lat[x],path.lon[x],4);
+
+ lat0=lat1;
+ lon0=lon1;
+
+ fgets(string,78,fd);
+
+ } while (strncmp(string,"END",3)!=0 && feof(fd)==0);
+
+ fgets(string,78,fd);
+
+ } while (strncmp(string,"END",3)!=0 && feof(fd)==0);
+
+ fclose(fd);
+
+ fprintf(stdout,"Done!\n");
+ fflush(stdout);
+ }
+ else
+ fprintf(stderr,"*** ERROR: \"%s\": not found!\n",filename);
+}
+
+void ReadLRParm(char *txsite_filename)
+{
+ /* This function reads Longley-Rice parameter data for the
+ transmitter site. The file name is the same as the txsite,
+ except the filename extension is .lrp. If the needed file
+ is not found, then the file "splat.lrp" is read from the
+ current working directory. Failure to load this file will
+ result in the default parameters hard coded into this
+ being used and written to "splat.lrp". */
+
+ double din;
+ char filename[255], lookup[256], string[80];
+ int iin, ok=0, x;
+ FILE *fd=NULL, *outfile=NULL;
+
+ /* Default parameters in case things go bad */
+
+ LR.eps_dielect=15.0;
+ LR.sgm_conductivity=0.005;
+ LR.eno_ns_surfref=301.0;
+ LR.frq_mhz=300.0;
+ LR.radio_climate=5;
+ LR.pol=0;
+ LR.conf=0.50;
+ LR.rel=0.50;
+
+ /* Modify txsite filename to one with a .lrp extension. */
+
+ strncpy(filename,txsite_filename,255);
+
+ for (x=0; filename[x]!='.' && filename[x]!=0 && filename[x]!='\n' && x<249; x++);
+
+ filename[x]='.';
+ filename[x+1]='l';
+ filename[x+2]='r';
+ filename[x+3]='p';
+ filename[x+4]=0;
+
+ /* Small lookup table to parse file, pass
+ numeric data, and ignore comments. */
+
+ for (x=0; x<=255; x++)
+ lookup[x]=0;
+
+ /* Valid characters */
+
+ for (x=48; x<=57; x++)
+ lookup[x]=x;
+
+ lookup[32]=32;
+ lookup[46]='.';
+
+ fd=fopen(filename,"r");
+
+ if (fd==NULL)
+ {
+ /* Load default "splat.lrp" file */
+
+ strncpy(filename,"splat.lrp\0",10);
+ fd=fopen(filename,"r");
+ }
+
+ if (fd!=NULL)
+ {
+ fgets(string,80,fd);
+
+ for (x=0; lookup[(int)string[x]] && x<20; x++)
+ string[x]=lookup[(int)string[x]];
+
+ string[x]=0;
+
+ ok=sscanf(string,"%lf", &din);
+
+ if (ok)
+ {
+ LR.eps_dielect=din;
+
+ fgets(string,80,fd);
+
+ for (x=0; lookup[(int)string[x]] && x<20; x++)
+ string[x]=lookup[(int)string[x]];
+
+ string[x]=0;
+
+ ok=sscanf(string,"%lf", &din);
+ }
+
+ if (ok)
+ {
+ LR.sgm_conductivity=din;
+
+ fgets(string,80,fd);
+
+ for (x=0; lookup[(int)string[x]] && x<20; x++)
+ string[x]=lookup[(int)string[x]];
+
+ string[x]=0;
+
+ ok=sscanf(string,"%lf", &din);
+ }
+
+ if (ok)
+ {
+ LR.eno_ns_surfref=din;
+
+ fgets(string,80,fd);
+
+ for (x=0; lookup[(int)string[x]] && x<20; x++)
+ string[x]=lookup[(int)string[x]];
+
+ string[x]=0;
+
+ ok=sscanf(string,"%lf", &din);
+ }
+
+ if (ok)
+ {
+ LR.frq_mhz=din;
+
+ fgets(string,80,fd);
+
+ for (x=0; lookup[(int)string[x]] && x<20; x++)
+ string[x]=lookup[(int)string[x]];
+
+ string[x]=0;
+
+ ok=sscanf(string,"%d", &iin);
+ }
+
+ if (ok)
+ {
+ LR.radio_climate=iin;
+
+ fgets(string,80,fd);
+
+ for (x=0; lookup[(int)string[x]] && x<20; x++)
+ string[x]=lookup[(int)string[x]];
+
+ string[x]=0;
+
+ ok=sscanf(string,"%d", &iin);
+ }
+
+ if (ok)
+ {
+ LR.pol=iin;
+
+ fgets(string,80,fd);
+
+ for (x=0; lookup[(int)string[x]] && x<20; x++)
+ string[x]=lookup[(int)string[x]];
+
+ string[x]=0;
+
+ ok=sscanf(string,"%lf", &din);
+ }
+
+ if (ok)
+ {
+ LR.conf=din;
+
+ fgets(string,80,fd);
+
+ for (x=0; lookup[(int)string[x]] && x<20; x++)
+ string[x]=lookup[(int)string[x]];
+
+ string[x]=0;
+
+ ok=sscanf(string,"%lf", &din);
+ }
+
+ if (ok)
+ LR.rel=din;
+
+ fclose(fd);
+ }
+
+ if (fd==NULL)
+ {
+ /* Create a "splat.lrp" file since one
+ could not be successfully loaded. */
+
+ outfile=fopen("splat.lrp","w");
+
+ fprintf(outfile,"%.3f\t; Earth Dielectric Constant (Relative permittivity)\n",LR.eps_dielect);
+
+ fprintf(outfile,"%.3f\t; Earth Conductivity (Siemens per meter)\n", LR.sgm_conductivity);
+
+ fprintf(outfile,"%.3f\t; Atmospheric Bending Constant (N-Units)\n",LR.eno_ns_surfref);
+
+ fprintf(outfile,"%.3f\t; Frequency in MHz (20 MHz to 20 GHz)\n", LR.frq_mhz);
+
+ fprintf(outfile,"%d\t; Radio Climate\n",LR.radio_climate);
+
+ fprintf(outfile,"%d\t; Polarization (0 = Horizontal, 1 = Vertical)\n", LR.pol);
+
+ fprintf(outfile,"%.2f\t; Fraction of situations\n",LR.conf);
+
+ fprintf(outfile, "%.2f\t; Fraction of time\n",LR.rel);
+
+ fprintf(outfile,"\nPlease consult SPLAT! documentation for the meaning and use of this data.\n");
+
+ fclose(outfile);
+
+ fprintf(stderr,"\n%c*** There were problems reading your \"%s\" file! ***\nA \"splat.lrp\" file was written to your directory with default data.\n",7,filename);
+ }
+
+ if (fd==NULL || ok==0)
+ fprintf(stderr,"Longley-Rice default parameters have been assumed for this analysis.\n");
+}
+
+struct site los(struct site source, struct site destination)
+{
+ /* This function determines whether a line-of-sight path
+ unobstructed by terrain exists between source (transmitter)
+ and destination (receiver) based on the geographical
+ locations of the two sites, their respective antenna
+ heights above ground, and the terrain between them.
+ A site structure is returned upon completion. If the
+ first character of site.name is ' ', then a clear path
+ exists between source and destination. If the first
+ character is '*', then an obstruction exists, and the
+ site.lat and site.lon elements of the structure provide
+ the geographical location of the obstruction. */
+
+ int x;
+ char block;
+ struct site test, blockage;
+ register double distance, tx_alt, rx_alt,
+ cos_xmtr_angle, cos_test_angle, test_alt;
+
+ ReadPath(source,destination);
+
+ distance=5280.0*Distance(source,destination);
+ tx_alt=earthradius+source.alt+GetElevation(source);
+ rx_alt=earthradius+destination.alt+GetElevation(destination);
+
+ /* Elevation angle of the xmtr (source) as seen by the rcvr */
+
+ cos_xmtr_angle=((rx_alt*rx_alt)+(distance*distance)-(tx_alt*tx_alt))/(2.0*rx_alt*distance);
+
+ /* Determine the elevation angle of each discrete location
+ along the path between the receiver and transmitter.
+
+ Since obstructions are more likely due to terrain effects
+ closest to the receiver rather than farther away, we start
+ looking for potential obstructions from the receiver's
+ location, and work our way towards the transmitter.
+ This loop is broken when the first obstruction is
+ detected. If we can travel all the way to the transmitter
+ without detecting an obstruction, then we have a clear
+ unobstructed path between transmitter and receiver. */
+
+ for (x=path.length-1, block=0; x>0 && block==0; x--)
+ {
+ /* Build a structure for each test
+ point along the path to be surveyed. */
+
+ test.lat=path.lat[x];
+ test.lon=path.lon[x];
+
+ /* Measure the distance between the
+ test point and the receiver locations */
+
+ distance=5280.0*Distance(test,destination);
+ test_alt=earthradius+path.elevation[x];
+
+ /* Determine the cosine of the elevation of the test
+ point as seen from the location of the receiver */
+
+ cos_test_angle=((rx_alt*rx_alt)+(distance*distance)-(test_alt*test_alt))/(2.0*rx_alt*distance);
+
+ /* If the elevation angle to the test point (as seen from
+ the receiver) is greater than the elevation angle to the
+ transmitter (as seen by the receiver), then we have a
+ path obstructed by terrain. Note: Since we're comparing
+ the cosines of these angles rather than the angles
+ themselves (eliminating the call to acos() saves
+ considerable time), the following "if" statement is
+ reversed from what it would normally be if the angles
+ were compared. */
+
+ if (cos_xmtr_angle>cos_test_angle)
+ {
+ block=1;
+ blockage.lat=path.lat[x];
+ blockage.lon=path.lon[x];
+ blockage.alt=path.elevation[x];
+ blockage.name[0]='*';
+ }
+ }
+
+ if (block==0)
+ {
+ blockage.lat=0.0;
+ blockage.lon=0.0;
+ blockage.alt=0.0;
+ blockage.name[0]=' ';
+ }
+
+ return blockage;
+}
+
+void PlotPath(struct site source, struct site destination, char mask_value)
+{
+ /* This function analyzes the path between the source and
+ destination locations. It determines which points along
+ the path have line-of-sight visibility to the source.
+ Points along with path having line-of-sight visibility
+ to the source at an AGL altitude equal to that of the
+ destination location are stored by setting bit 1 in the
+ mask[][] array, which are displayed in green when PPM
+ maps are later generated by SPLAT!. */
+
+ char block;
+ int x, y;
+ register double cos_xmtr_angle, cos_test_angle, test_alt;
+ double distance, rx_alt, tx_alt;
+
+ ReadPath(source,destination);
+
+ for (y=0; y<path.length; y++)
+ {
+ /* Test this point only if it hasn't been already
+ tested and found to be free of obstructions. */
+
+ if ((GetMask(path.lat[y],path.lon[y])&mask_value)==0)
+ {
+ distance=5280.0*path.distance[y];
+ tx_alt=earthradius+source.alt+path.elevation[0];
+ rx_alt=earthradius+destination.alt+path.elevation[y];
+
+ /* Calculate the cosine of the elevation of the
+ transmitter as seen at the temp rx point. */
+
+ cos_xmtr_angle=((rx_alt*rx_alt)+(distance*distance)-(tx_alt*tx_alt))/(2.0*rx_alt*distance);
+
+ for (x=y, block=0; x>=0 && block==0; x--)
+ {
+ distance=5280.0*(path.distance[y]-path.distance[x]);
+ test_alt=earthradius+path.elevation[x];
+
+ cos_test_angle=((rx_alt*rx_alt)+(distance*distance)-(test_alt*test_alt))/(2.0*rx_alt*distance);
+
+ /* Compare these two angles to determine if
+ a blockage exists. Since we're comparing
+ the cosines of these angles rather than
+ the angles themselves, the following "if"
+ statement is reversed from what it would
+ be if the actual angles were compared. */
+
+ if (cos_xmtr_angle>cos_test_angle)
+ block=1;
+ }
+
+ if (block==0)
+ OrMask(path.lat[y],path.lon[y],mask_value);
+ }
+ }
+}
+
+void PlotLRPath(struct site source, struct site destination)
+{
+ /* This function plots the RF signal path loss
+ between source and destination points based
+ on the Longley-Rice propagation model. */
+
+ char strmode[100];
+ int x, y, errnum;
+ double loss;
+
+ ReadPath(source,destination);
+ elev_l[1]=0.04*METERS_PER_MILE;
+
+ for (x=0; x<path.length; x++)
+ elev_l[x+2]=path.elevation[x]*METERS_PER_FOOT;
+
+ for (y=0; y<path.length; y++)
+ {
+ /* Test this point only if it
+ has not already been tested. */
+
+ if (GetMask(path.lat[y],path.lon[y])==0 && 0.04*y<=max_range)
+ {
+ elev_l[0]=y+1;
+
+ point_to_point(elev_l,source.alt*METERS_PER_FOOT,
+ destination.alt*METERS_PER_FOOT,
+ LR.eps_dielect, LR.sgm_conductivity, LR.eno_ns_surfref,
+ LR.frq_mhz, LR.radio_climate, LR.pol, LR.conf, LR.rel,
+ loss, strmode, errnum);
+
+ /* Note: PASS BY REFERENCE ... loss and errnum are pass
+ by reference, only used in this file by this function */
+
+ if (loss>225.0)
+ loss=225.0;
+
+ if (loss<75.0)
+ loss=75.0;
+
+ loss-=75.0;
+ loss/=10.0;
+ loss+=1.0;
+
+ OrMask(path.lat[y],path.lon[y],((unsigned char)(loss))<<3);
+ }
+
+ else if (GetMask(path.lat[y],path.lon[y])==0 && 0.04*y>max_range)
+ OrMask(path.lat[y],path.lon[y],1);
+ }
+}
+
+void PlotCoverage(struct site source, double altitude)
+{
+ /* This function performs a 360 degree sweep around the
+ transmitter site (source location), and plots the
+ line-of-sight coverage of the transmitter on the SPLAT!
+ generated topographic map based on a receiver located
+ at the specified altitude (in feet AGL). Results
+ are stored in memory, and written out in the form
+ of a topographic map when the WritePPM() function
+ is later invoked. */
+
+ double lat, lon, one_pixel;
+ static unsigned char mask_value;
+ int z, count;
+ struct site edge;
+ unsigned char symbol[4], x;
+
+ /* Initialize mask_value */
+
+ if (mask_value!=8 && mask_value!=16 && mask_value!=32)
+ mask_value=1;
+
+ one_pixel=1.0/1200.0;
+
+ symbol[0]='.';
+ symbol[1]='o';
+ symbol[2]='O';
+ symbol[3]='o';
+
+ count=0;
+
+ fprintf(stdout,"\nComputing line-of-sight coverage of %s with an RX antenna\nat %.2f feet AGL:\n\n 0%c to 25%c ",source.name,altitude,37,37);
+ fflush(stdout);
+
+ /* 18.75=1200 pixels/degree divided by 64 loops
+ per progress indicator symbol (.oOo) printed. */
+
+ z=(int)(18.75*(max_west-min_west));
+
+ for (lon=min_west, x=0; lon<=max_west; lon+=one_pixel)
+ {
+ edge.lat=max_north;
+ edge.lon=lon;
+ edge.alt=altitude;
+
+ PlotPath(source,edge,mask_value);
+ count++;
+
+ if (count==z)
+ {
+ fprintf(stdout,"%c",symbol[x]);
+ fflush(stdout);
+ count=0;
+
+ if (x==3)
+ x=0;
+ else
+ x++;
+ }
+ }
+
+ count=0;
+ fprintf(stdout,"\n25%c to 50%c ",37,37);
+ fflush(stdout);
+
+ z=(int)(18.75*(max_north-min_north));
+
+ for (lat=max_north, x=0; lat>=min_north; lat-=one_pixel)
+ {
+ edge.lat=lat;
+ edge.lon=min_west;
+ edge.alt=altitude;
+
+ PlotPath(source,edge,mask_value);
+ count++;
+
+ if (count==z)
+ {
+ fprintf(stdout,"%c",symbol[x]);
+ fflush(stdout);
+ count=0;
+
+ if (x==3)
+ x=0;
+ else
+ x++;
+ }
+ }
+
+ count=0;
+ fprintf(stdout,"\n50%c to 75%c ",37,37);
+ fflush(stdout);
+
+ z=(int)(18.75*(max_west-min_west));
+
+ for (lon=min_west, x=0; lon<=max_west; lon+=one_pixel)
+ {
+ edge.lat=min_north;
+ edge.lon=lon;
+ edge.alt=altitude;
+
+ PlotPath(source,edge,mask_value);
+ count++;
+
+ if (count==z)
+ {
+ fprintf(stdout,"%c",symbol[x]);
+ fflush(stdout);
+ count=0;
+
+ if (x==3)
+ x=0;
+ else
+ x++;
+ }
+ }
+
+ count=0;
+ fprintf(stdout,"\n75%c to 100%c ",37,37);
+ fflush(stdout);
+
+ z=(int)(18.75*(max_north-min_north));
+
+ for (lat=min_north, x=0; lat<=max_north; lat+=one_pixel)
+ {
+ edge.lat=lat;
+ edge.lon=max_west;
+ edge.alt=altitude;
+
+ PlotPath(source,edge,mask_value);
+ count++;
+
+ if (count==z)
+ {
+ fprintf(stdout,"%c",symbol[x]);
+ fflush(stdout);
+ count=0;
+
+ if (x==3)
+ x=0;
+ else
+ x++;
+ }
+ }
+
+ fprintf(stdout,"\nDone!\n");
+ fflush(stdout);
+
+ /* Assign next mask value */
+
+ switch (mask_value)
+ {
+ case 1:
+ mask_value=8;
+ break;
+
+ case 8:
+ mask_value=16;
+ break;
+
+ case 16:
+ mask_value=32;
+ }
+}
+
+void PlotLRMap(struct site source, double altitude)
+{
+ /* This function performs a 360 degree sweep around the
+ transmitter site (source location), and plots the
+ Longley-Rice attenuation on the SPLAT! generated
+ topographic map based on a receiver located at
+ the specified altitude (in feet AGL). Results
+ are stored in memory, and written out in the form
+ of a topographic map when the WritePPMLR() function
+ is later invoked. */
+
+ int z, count;
+ struct site edge;
+ double lat, lon, one_pixel;
+ unsigned char symbol[4], x;
+
+ one_pixel=1.0/1200.0;
+
+ symbol[0]='.';
+ symbol[1]='o';
+ symbol[2]='O';
+ symbol[3]='o';
+
+ count=0;
+
+ fprintf(stdout,"\nComputing Longley-Rice coverage of %s ", source.name);
+ fprintf(stdout,"out to a radius\nof %.2f miles with an RX antenna at %.2f feet AGL:\n\n 0%c to 25%c ",max_range,altitude,37,37);
+ fflush(stdout);
+
+ /* 18.75=1200 pixels/degree divided by 64 loops
+ per progress indicator symbol (.oOo) printed. */
+
+ z=(int)(18.75*(max_west-min_west));
+
+ for (lon=min_west, x=0; lon<=max_west; lon+=one_pixel)
+ {
+ edge.lat=max_north;
+ edge.lon=lon;
+ edge.alt=altitude;
+
+ PlotLRPath(source,edge);
+ count++;
+
+ if (count==z)
+ {
+ fprintf(stdout,"%c",symbol[x]);
+ fflush(stdout);
+ count=0;
+
+ if (x==3)
+ x=0;
+ else
+ x++;
+ }
+ }
+
+ count=0;
+ fprintf(stdout,"\n25%c to 50%c ",37,37);
+ fflush(stdout);
+
+ z=(int)(18.75*(max_north-min_north));
+
+ for (lat=max_north, x=0; lat>=min_north; lat-=one_pixel)
+ {
+ edge.lat=lat;
+ edge.lon=min_west;
+ edge.alt=altitude;
+
+ PlotLRPath(source,edge);
+ count++;
+
+ if (count==z)
+ {
+ fprintf(stdout,"%c",symbol[x]);
+ fflush(stdout);
+ count=0;
+
+ if (x==3)
+ x=0;
+ else
+ x++;
+ }
+ }
+
+ count=0;
+ fprintf(stdout,"\n50%c to 75%c ",37,37);
+ fflush(stdout);
+
+ z=(int)(18.75*(max_west-min_west));
+
+ for (lon=min_west, x=0; lon<=max_west; lon+=one_pixel)
+ {
+ edge.lat=min_north;
+ edge.lon=lon;
+ edge.alt=altitude;
+
+ PlotLRPath(source,edge);
+ count++;
+
+ if (count==z)
+ {
+ fprintf(stdout,"%c",symbol[x]);
+ fflush(stdout);
+ count=0;
+
+ if (x==3)
+ x=0;
+ else
+ x++;
+ }
+ }
+
+ count=0;
+ fprintf(stdout,"\n75%c to 100%c ",37,37);
+ fflush(stdout);
+
+ z=(int)(18.75*(max_north-min_north));
+
+ for (lat=min_north, x=0; lat<=max_north; lat+=one_pixel)
+ {
+ edge.lat=lat;
+ edge.lon=max_west;
+ edge.alt=altitude;
+
+ PlotLRPath(source,edge);
+ count++;
+
+ if (count==z)
+ {
+ fprintf(stdout,"%c",symbol[x]);
+ fflush(stdout);
+ count=0;
+
+ if (x==3)
+ x=0;
+ else
+ x++;
+ }
+ }
+
+ fprintf(stdout,"\nDone!\n");
+ fflush(stdout);
+}
+
+void WritePPM(char *filename)
+{
+ /* This function generates a topographic map in Portable Pix Map
+ (PPM) format based on logarithmically scaled topology data,
+ as well as the content of flags held in the mask[][] array.
+ The image created is rotated counter-clockwise 90 degrees
+ from its representation in dem[][] so that north points
+ up and east points right in the image generated. */
+
+ int indx, x, x0, y0, minlat, minlon;
+ unsigned width, height, output;
+ unsigned char found, mask;
+ char mapfile[255];
+ double conversion, lat, lon, one_over_gamma, one_pixel;
+ FILE *fd;
+
+ one_pixel=1.0/1200.0;
+ one_over_gamma=1.0/GAMMA;
+ conversion=255.0/pow((double)(max_elevation-min_elevation),one_over_gamma);
+
+ width=1200*(max_west-min_west);
+ height=1200*(max_north-min_north);
+
+ if (filename[0]==0)
+ strncpy(mapfile, "map.ppm\0",8);
+ else
+ {
+ for (x=0; filename[x]!='.' && filename[x]!=0 && x<250; x++)
+ mapfile[x]=filename[x];
+
+ mapfile[x]='.';
+ mapfile[x+1]='p';
+ mapfile[x+2]='p';
+ mapfile[x+3]='m';
+ mapfile[x+4]=0;
+ }
+
+ fd=fopen(mapfile,"wb");
+
+ fprintf(fd,"P6\n%u %u\n255\n",width,height);
+
+ fprintf(stdout,"\nWriting \"%s\" (%ux%u pixmap image)... ",mapfile,width,height);
+ fflush(stdout);
+
+ for (lat=(double)max_north; lat>=(double)min_north; lat-=one_pixel)
+ {
+ for (lon=(double)max_west; lon>=(double)min_west; lon-=one_pixel)
+ {
+ minlat=(int)floor(lat);
+ minlon=(int)floor(lon);
+
+ for (indx=0, found=0; indx<MAXSLOTS && found==0;)
+ if (minlat==dem[indx].min_north && minlon==dem[indx].min_west)
+ found=1;
+ else
+ indx++;
+ if (found)
+ {
+ x0=(int)((1199.0*(lat-floor(lat)))+0.5);
+ y0=(int)((1199.0*(lon-floor(lon)))+0.5);
+
+ mask=dem[indx].mask[x0][y0];
+
+ if (mask&2)
+ /* Text Labels: Red */
+ fprintf(fd,"%c%c%c",255,0,0);
+
+ else if (mask&4)
+ /* County Boundaries: Light Cyan */
+ fprintf(fd,"%c%c%c",128,128,255);
+
+ else switch (mask&57)
+ {
+ case 1:
+ /* TX1: Green */
+ fprintf(fd,"%c%c%c",0,255,0);
+ break;
+
+ case 8:
+ /* TX2: Cyan */
+ fprintf(fd,"%c%c%c",0,255,255);
+ break;
+
+ case 9:
+ /* TX1 + TX2: Yellow */
+ fprintf(fd,"%c%c%c",255,255,0);
+ break;
+
+ case 16:
+ /* TX3: Medium Violet */
+ fprintf(fd,"%c%c%c",147,112,219);
+ break;
+
+ case 17:
+ /* TX1 + TX3: Pink */
+ fprintf(fd,"%c%c%c",255,192,203);
+ break;
+
+ case 24:
+ /* TX2 + TX3: Orange */
+ fprintf(fd,"%c%c%c",255,165,0);
+ break;
+
+ case 25:
+ /* TX1 + TX2 + TX3: Dark Green */
+ fprintf(fd,"%c%c%c",0,100,0);
+ break;
+
+ case 32:
+ /* TX4: Sienna 1 */
+ fprintf(fd,"%c%c%c",255,130,71);
+ break;
+
+ case 33:
+ /* TX1 + TX4: Green Yellow */
+ fprintf(fd,"%c%c%c",173,255,47);
+ break;
+
+ case 40:
+ /* TX2 + TX4: Dark Sea Green 1 */
+ fprintf(fd,"%c%c%c",193,255,193);
+ break;
+
+ case 41:
+ /* TX1 + TX2 + TX4: Blanched Almond */
+ fprintf(fd,"%c%c%c",255,235,205);
+ break;
+
+ case 48:
+ /* TX3 + TX4: Dark Turquoise */
+ fprintf(fd,"%c%c%c",0,206,209);
+ break;
+
+ case 49:
+ /* TX1 + TX3 + TX4: Medium Spring Green */
+ fprintf(fd,"%c%c%c",0,250,154);
+ break;
+
+ case 56:
+ /* TX2 + TX3 + TX4: Tan */
+ fprintf(fd,"%c%c%c",210,180,140);
+ break;
+
+ case 57:
+ /* TX1 + TX2 + TX3 + TX4: Gold2 */
+ fprintf(fd,"%c%c%c",238,201,0);
+ break;
+
+ default:
+ /* Water: Medium Blue */
+ if (dem[indx].data[x0][y0]==0)
+ fprintf(fd,"%c%c%c",0,0,170);
+ else
+ {
+ /* Elevation: Greyscale */
+ output=(unsigned)(0.5+pow((double)(dem[indx].data[x0][y0]-min_elevation),one_over_gamma)*conversion);
+ fprintf(fd,"%c%c%c",output,output,output);
+ }
+ }
+ }
+ }
+ }
+
+ fclose(fd);
+ fprintf(stdout,"Done!\n");
+ fflush(stdout);
+}
+
+void WritePPMLR(char *filename)
+{
+ /* This function generates a topographic map in Portable Pix Map
+ (PPM) format based on the content of flags held in the mask[][]
+ array (only). The image created is rotated counter-clockwise
+ 90 degrees from its representation in dem[][] so that north
+ points up and east points right in the image generated. */
+
+ int indx, x, t, t2, x0, y0, minlat, minlon;
+ unsigned width, height, output;
+ unsigned char found, mask;
+ char mapfile[255];
+ double conversion, lat, lon, one_over_gamma, one_pixel;
+ FILE *fd;
+
+ one_pixel=1.0/1200.0;
+ one_over_gamma=1.0/GAMMA;
+ conversion=255.0/pow((double)(max_elevation-min_elevation),one_over_gamma);
+
+ width=1200*(max_west-min_west);
+ height=1200*(max_north-min_north);
+
+ if (filename[0]==0)
+ strncpy(mapfile, "map.ppm\0",8);
+ else
+ {
+ for (x=0; filename[x]!='.' && filename[x]!=0 && x<250; x++)
+ mapfile[x]=filename[x];
+
+ mapfile[x]='.';
+ mapfile[x+1]='p';
+ mapfile[x+2]='p';
+ mapfile[x+3]='m';
+ mapfile[x+4]=0;
+ }
+
+ fd=fopen(mapfile,"wb");
+
+ fprintf(fd,"P6\n%u %u\n255\n",width,height+30);
+
+ fprintf(stdout,"\nWriting \"%s\" (%ux%u pixmap image)... ",mapfile,width,height+30);
+ fflush(stdout);
+
+ for (lat=(double)max_north; lat>=(double)min_north; lat-=one_pixel)
+ {
+ for (lon=(double)max_west; lon>=(double)min_west; lon-=one_pixel)
+ {
+ minlat=(int)floor(lat);
+ minlon=(int)floor(lon);
+
+ for (indx=0, found=0; indx<MAXSLOTS && found==0;)
+ if (minlat==dem[indx].min_north && minlon==dem[indx].min_west)
+ found=1;
+ else
+ indx++;
+ if (found)
+ {
+ x0=(int)((1199.0*(lat-floor(lat)))+0.5);
+ y0=(int)((1199.0*(lon-floor(lon)))+0.5);
+
+ mask=dem[indx].mask[x0][y0];
+
+ if (mask&2)
+ {
+ /* Text Labels - Black or Red */
+ if (mask&120)
+ fprintf(fd,"%c%c%c",0,0,0);
+ else
+ fprintf(fd,"%c%c%c",255,0,0);
+ }
+
+ else if (mask&4)
+ /* County Boundaries: Black */
+ fprintf(fd,"%c%c%c",0,0,0);
+
+ else if (mask&1 && !((mask&248)==192))
+ {
+ /* Outside Analysis Range */
+ /* Display Greyscale / Sea Level */
+
+ if (dem[indx].data[x0][y0]==0)
+ fprintf(fd,"%c%c%c",0,0,170);
+ else
+ {
+ output=(unsigned)(0.5+pow((double)(dem[indx].data[x0][y0]-min_elevation),one_over_gamma)*conversion);
+ fprintf(fd,"%c%c%c",output,output,output);
+ }
+ }
+
+ else switch ((mask&248)>>3)
+ {
+ case 0:
+ /* Inside range, but no coverage.
+ Display Sea Level / Terrain */
+
+ if (dem[indx].data[x0][y0]==0)
+ fprintf(fd,"%c%c%c",0,0,170);
+ else
+ {
+ /* Elevation: Greyscale */
+ output=(unsigned)(0.5+pow((double)(dem[indx].data[x0][y0]-min_elevation),one_over_gamma)*conversion);
+ fprintf(fd,"%c%c%c",output,output,output);
+ }
+
+ break;
+
+ case 1:
+ /* Green */
+ fprintf(fd,"%c%c%c",0,255,0);
+ break;
+
+ case 2:
+ /* Pink */
+ fprintf(fd,"%c%c%c",255,192,203);
+ break;
+
+ case 3:
+ /* Cyan */
+ fprintf(fd,"%c%c%c",0,255,255);
+ break;
+
+ case 4:
+ /* Yellow */
+ fprintf(fd,"%c%c%c",255,255,0);
+ break;
+
+ case 5:
+ /* Medium Violet */
+ fprintf(fd,"%c%c%c",161,131,224);
+ break;
+
+ case 6:
+ /* Orange */
+ fprintf(fd,"%c%c%c",255,165,0);
+ break;
+
+ case 7:
+ /* Light Green */
+ fprintf(fd,"%c%c%c",193,255,193);
+ break;
+
+ case 8:
+ /* Red Pink */
+ fprintf(fd,"%c%c%c",255,108,108);
+ break;
+
+ case 9:
+ /* TX1 + TX4: Green Yellow */
+ fprintf(fd,"%c%c%c",173,255,47);
+ break;
+
+ case 10:
+ /* Blanched Almond */
+ fprintf(fd,"%c%c%c",255,235,184);
+ break;
+
+ case 11:
+ /* Dark Turquoise */
+ fprintf(fd,"%c%c%c",0,206,209);
+ break;
+
+ case 12:
+ /* Tan */
+ fprintf(fd,"%c%c%c",210,180,140);
+ break;
+
+ case 13:
+ /* Magenta 1 */
+ fprintf(fd,"%c%c%c",243,110,205);
+ break;
+
+ case 14:
+ /* Gold2 */
+ fprintf(fd,"%c%c%c",238,201,0);
+ break;
+
+ case 15:
+ /* Medium Spring Green */
+ fprintf(fd,"%c%c%c",0,250,154);
+ break;
+
+ case 16:
+ /* Very light Blue */
+ fprintf(fd,"%c%c%c",244,244,255);
+ break;
+
+ default:
+ /* Land / Sea */
+ if (dem[indx].data[x0][y0]==0)
+ fprintf(fd,"%c%c%c",0,0,170);
+ else
+ {
+ /* Elevation: Greyscale */
+ output=(unsigned)(0.5+pow((double)(dem[indx].data[x0][y0]-min_elevation),one_over_gamma)*conversion);
+ fprintf(fd,"%c%c%c",output,output,output);
+ }
+ }
+ }
+ }
+ }
+
+ x0=width/16;
+
+ for (y0=0; y0<30; y0++)
+ {
+ for (indx=0; indx<16; indx++)
+ {
+ for (x=0; x<x0; x++)
+ {
+ t=indx;
+ t2=indx+8;
+
+ if (y0>=10 && y0<=18)
+ {
+ if (t2>9)
+ {
+ if (x>=11 && x<=17)
+ if (smallfont[t2/10][y0-10][x-11])
+ t=255;
+ }
+
+ if (x>=19 && x<=25)
+ if (smallfont[t2%10][y0-10][x-19])
+ t=255;
+
+ if (x>=27 && x<=33)
+ if (smallfont[0][y0-10][x-27])
+ t=255;
+ }
+
+ switch (t)
+ {
+ case 0:
+ /* Green */
+ fprintf(fd,"%c%c%c",0,255,0);
+ break;
+
+ case 1:
+ /* Pink */
+ fprintf(fd,"%c%c%c",255,192,203);
+ break;
+
+ case 2:
+ /* Cyan */
+ fprintf(fd,"%c%c%c",0,255,255);
+ break;
+
+ case 3:
+ /* Yellow */
+ fprintf(fd,"%c%c%c",255,255,0);
+ break;
+
+ case 4:
+ /* Medium Violet */
+ fprintf(fd,"%c%c%c",161,131,224);
+ break;
+
+ case 5:
+ /* Orange */
+ fprintf(fd,"%c%c%c",255,165,0);
+ break;
+
+ case 6:
+ /* Light Green */
+ fprintf(fd,"%c%c%c",193,255,193);
+ break;
+
+ case 7:
+ /* Red Pink */
+ fprintf(fd,"%c%c%c",255,108,108);
+ break;
+
+ case 8:
+ /* Green Yellow */
+ fprintf(fd,"%c%c%c",173,255,47);
+ break;
+
+ case 9:
+ /* Blanched Almond */
+ fprintf(fd,"%c%c%c",255,235,184);
+ break;
+
+ case 10:
+ /* Dark Turquoise */
+ fprintf(fd,"%c%c%c",0,206,209);
+ break;
+
+ case 11:
+ /* Tan */
+ fprintf(fd,"%c%c%c",210,180,140);
+ break;
+
+ case 12:
+ /* Magenta 1 */
+ fprintf(fd,"%c%c%c",243,110,205);
+ break;
+
+ case 13:
+ /* Gold2 */
+ fprintf(fd,"%c%c%c",238,201,0);
+ break;
+
+ case 14:
+ /* Medium Spring Green */
+ fprintf(fd,"%c%c%c",0,250,154);
+ break;
+
+ case 255:
+ /* Black */
+ fprintf(fd,"%c%c%c",0,0,0);
+ break;
+
+ default:
+ /* Very Light Blue */
+ fprintf(fd,"%c%c%c",244,244,255);
+ }
+ }
+ }
+ }
+
+ fclose(fd);
+ fprintf(stdout,"Done!\n");
+ fflush(stdout);
+}
+
+void GraphTerrain(struct site source, struct site destination, char *name)
+{
+ /* This function invokes gnuplot to generate an appropriate
+ output file indicating the terrain profile between the source
+ and destination locations. "filename" is the name assigned
+ to the output file generated by gnuplot. The filename extension
+ is used to set gnuplot's terminal setting and output file type.
+ If no extension is found, .gif is assumed. */
+
+ int x, y, z;
+ char filename[255], term[15], ext[15];
+ FILE *fd=NULL;
+
+ ReadPath(destination,source);
+
+ fd=fopen("profile.gp","wb");
+
+ for (x=0; x<path.length; x++)
+ fprintf(fd,"%f\t%f\n",path.distance[x],path.elevation[x]);
+
+ fclose(fd);
+
+ if (name[0]==0)
+ {
+ /* Default filename and output file type */
+
+ strncpy(filename,"profile\0",8);
+ strncpy(term,"gif\0",4);
+ strncpy(ext,"gif\0",4);
+ }
+
+ else
+ {
+ /* Grab extension and terminal type from "name" */
+
+ for (x=0; name[x]!='.' && name[x]!=0 && x<254; x++)
+ filename[x]=name[x];
+
+ if (name[x]=='.')
+ {
+ for (y=0, z=x, x++; name[x]!=0 && x<254 && y<14; x++, y++)
+ {
+ term[y]=tolower(name[x]);
+ ext[y]=term[y];
+ }
+
+ ext[y]=0;
+ term[y]=0;
+ filename[z]=0;
+ }
+
+ else
+ { /* No extension -- Default is gif */
+
+ filename[x]=0;
+ strncpy(term,"gif\0",4);
+ strncpy(ext,"gif\0",4);
+ }
+ }
+
+ /* Either .ps or .postscript may be used
+ as an extension for postscript output. */
+
+ if (strncmp(term,"postscript",10)==0)
+ strncpy(ext,"ps\0",3);
+
+ else if (strncmp(ext,"ps",2)==0)
+ strncpy(term,"postscript\0",11);
+
+ fprintf(stdout,"Writing \"%s.%s\"...",filename,ext);
+ fflush(stdout);
+
+ fd=fopen("splat.gp","w");
+ fprintf(fd,"set grid\n");
+ fprintf(fd,"set autoscale\n");
+ fprintf(fd,"set term %s\n",term);
+ fprintf(fd,"set title \"SPLAT! Terrain Profile\"\n");
+ fprintf(fd,"set xlabel \"Distance Between %s and %s (miles)\"\n",destination.name,source.name);
+ fprintf(fd,"set ylabel \"Ground Elevation Above Sea Level (feet)\"\n");
+ fprintf(fd,"set output \"%s.%s\"\n",filename,ext);
+ fprintf(fd,"plot \"profile.gp\" title \"\" with lines\n");
+ fclose(fd);
+
+ x=system("gnuplot splat.gp");
+
+ if (x!=-1)
+ {
+ unlink("splat.gp");
+ unlink("profile.gp");
+ fprintf(stdout," Done!\n");
+ fflush(stdout);
+ }
+
+ else
+ fprintf(stderr,"\n*** ERROR: Error occurred invoking gnuplot!\n");
+}
+
+void GraphElevation(struct site source, struct site destination, char *name)
+{
+ /* This function invokes gnuplot to generate an appropriate
+ output file indicating the terrain profile between the source
+ and destination locations. "filename" is the name assigned
+ to the output file generated by gnuplot. The filename extension
+ is used to set gnuplot's terminal setting and output file type.
+ If no extension is found, .gif is assumed. */
+
+ int x, y, z;
+ char filename[255], term[15], ext[15];
+ double angle, refangle, maxangle=-90.0;
+ struct site remote;
+ FILE *fd=NULL, *fd2=NULL;
+
+ ReadPath(destination,source); /* destination=RX, source=TX */
+ refangle=ElevationAngle(destination,source);
+
+ fd=fopen("profile.gp","wb");
+ fd2=fopen("reference.gp","wb");
+
+ for (x=1; x<path.length-1; x++)
+ {
+ remote.lat=path.lat[x];
+ remote.lon=path.lon[x];
+ remote.alt=0.0;
+ angle=ElevationAngle(destination,remote);
+ fprintf(fd,"%f\t%f\n",path.distance[x],angle);
+ fprintf(fd2,"%f\t%f\n",path.distance[x],refangle);
+
+ if (angle>maxangle)
+ maxangle=angle;
+ }
+
+ fprintf(fd,"%f\t%f\n",path.distance[path.length-1],refangle);
+ fprintf(fd2,"%f\t%f\n",path.distance[path.length-1],refangle);
+
+ fclose(fd);
+ fclose(fd2);
+
+ if (name[0]==0)
+ {
+ /* Default filename and output file type */
+
+ strncpy(filename,"profile\0",8);
+ strncpy(term,"gif\0",4);
+ strncpy(ext,"gif\0",4);
+ }
+
+ else
+ {
+ /* Grab extension and terminal type from "name" */
+
+ for (x=0; name[x]!='.' && name[x]!=0 && x<254; x++)
+ filename[x]=name[x];
+
+ if (name[x]=='.')
+ {
+ for (y=0, z=x, x++; name[x]!=0 && x<254 && y<14; x++, y++)
+ {
+ term[y]=tolower(name[x]);
+ ext[y]=term[y];
+ }
+
+ ext[y]=0;
+ term[y]=0;
+ filename[z]=0;
+ }
+
+ else
+ { /* No extension -- Default is gif */
+
+ filename[x]=0;
+ strncpy(term,"gif\0",4);
+ strncpy(ext,"gif\0",4);
+ }
+ }
+
+ /* Either .ps or .postscript may be used
+ as an extension for postscript output. */
+
+ if (strncmp(term,"postscript",10)==0)
+ strncpy(ext,"ps\0",3);
+
+ else if (strncmp(ext,"ps",2)==0)
+ strncpy(term,"postscript\0",11);
+
+ fprintf(stdout,"Writing \"%s.%s\"...",filename,ext);
+ fflush(stdout);
+
+ fd=fopen("splat.gp","w");
+
+ fprintf(fd,"set grid\n");
+ fprintf(fd,"set yrange [%2.3f to %2.3f]\n", (-fabs(refangle)-0.25), maxangle+0.25);
+ fprintf(fd,"set term %s\n",term);
+ fprintf(fd,"set title \"SPLAT! Elevation Profile\"\n");
+ fprintf(fd,"set xlabel \"Distance Between %s and %s (miles)\"\n",destination.name,source.name);
+ fprintf(fd,"set ylabel \"Elevation Angle Along Path Between %s and %s (degrees)\"\n",destination.name,source.name);
+ fprintf(fd,"set output \"%s.%s\"\n",filename,ext);
+ fprintf(fd,"plot \"profile.gp\" title \"Real Earth Profile\" with lines, \"reference.gp\" title \"Line Of Sight Path\" with lines\n");
+
+ fclose(fd);
+
+ x=system("gnuplot splat.gp");
+
+ if (x!=-1)
+ {
+ unlink("splat.gp");
+ unlink("profile.gp");
+ unlink("reference.gp");
+
+ fprintf(stdout," Done!\n");
+ fflush(stdout);
+ }
+
+ else
+ fprintf(stderr,"\n*** ERROR: Error occurred invoking gnuplot!\n");
+}
+
+void GraphHeight(struct site source, struct site destination, char *name)
+{
+ /* This function invokes gnuplot to generate an appropriate
+ output file indicating the terrain profile between the source
+ and destination locations. What is plotted is the height of land
+ above or below a straight line between the receibe and transmit
+ sites. "filename" is the name assigned to the output file
+ generated by gnuplot. The filename extension is used to
+ set gnuplot's terminal setting and output file type.
+ If no extension is found, .gif is assumed. */
+
+ int x, y, z;
+ char filename[255], term[15], ext[15];
+ double a, b, c, height, refangle, cangle, maxheight=-100000.0,
+ minheight=100000.0;
+ struct site remote;
+ FILE *fd=NULL, *fd2=NULL;
+
+ ReadPath(destination,source); /* destination=RX, source=TX */
+ refangle=ElevationAngle(destination,source);
+ b=GetElevation(destination)+destination.alt+earthradius;
+
+ fd=fopen("profile.gp","wb");
+ fd2=fopen("reference.gp","wb");
+
+ for (x=1; x<path.length-1; x++)
+ {
+ remote.lat=path.lat[x];
+ remote.lon=path.lon[x];
+ remote.alt=0.0;
+
+ a=GetElevation(remote)+earthradius;
+
+ cangle=5280.0*Distance(destination,remote)/earthradius;
+
+ c=b*sin(refangle*deg2rad+HALFPI)/sin(HALFPI-refangle*deg2rad-cangle);
+
+ height=a-c;
+
+ fprintf(fd,"%f\t%f\n",path.distance[x],height);
+ fprintf(fd2,"%f\t%f\n",path.distance[x],0.);
+
+ if (height>maxheight)
+ maxheight=height;
+
+ if (height<minheight)
+ minheight=height;
+ }
+
+ fprintf(fd,"%f\t%f\n",path.distance[path.length-1],0.0);
+ fprintf(fd2,"%f\t%f\n",path.distance[path.length-1],0.0);
+
+ fclose(fd);
+ fclose(fd2);
+
+ if (name[0]==0)
+ {
+ /* Default filename and output file type */
+
+ strncpy(filename,"height\0",8);
+ strncpy(term,"gif\0",4);
+ strncpy(ext,"gif\0",4);
+ }
+
+ else
+ {
+ /* Grab extension and terminal type from "name" */
+
+ for (x=0; name[x]!='.' && name[x]!=0 && x<254; x++)
+ filename[x]=name[x];
+
+ if (name[x]=='.')
+ {
+ for (y=0, z=x, x++; name[x]!=0 && x<254 && y<14; x++, y++)
+ {
+ term[y]=tolower(name[x]);
+ ext[y]=term[y];
+ }
+
+ ext[y]=0;
+ term[y]=0;
+ filename[z]=0;
+ }
+
+ else
+ { /* No extension -- Default is gif */
+
+ filename[x]=0;
+ strncpy(term,"gif\0",4);
+ strncpy(ext,"gif\0",4);
+ }
+ }
+
+ /* Either .ps or .postscript may be used
+ as an extension for postscript output. */
+
+ if (strncmp(term,"postscript",10)==0)
+ strncpy(ext,"ps\0",3);
+
+ else if (strncmp(ext,"ps",2)==0)
+ strncpy(term,"postscript\0",11);
+
+ fprintf(stdout,"Writing \"%s.%s\"...",filename,ext);
+ fflush(stdout);
+
+ fd=fopen("splat.gp","w");
+
+ minheight-=20.0;
+ maxheight+=20.0;
+
+ if (maxheight<20.0)
+ maxheight=20.0;
+
+ fprintf(fd,"set grid\n");
+ fprintf(fd,"set yrange [%2.3f to %2.3f]\n", minheight, maxheight);
+ fprintf(fd,"set term %s\n",term);
+ fprintf(fd,"set title \"SPLAT! Height Profile\"\n");
+ fprintf(fd,"set xlabel \"Distance Between %s and %s (miles)\"\n",destination.name,source.name);
+ fprintf(fd,"set ylabel \"Ground Height Above Path Between %s and %s (feet)\"\n",destination.name,source.name);
+ fprintf(fd,"set output \"%s.%s\"\n",filename,ext);
+ fprintf(fd,"plot \"profile.gp\" title \"Real Earth Profile\" with lines, \"reference.gp\" title \"Line Of Sight Path\" with lines\n");
+
+ fclose(fd);
+
+ x=system("gnuplot splat.gp");
+
+ if (x!=-1)
+ {
+ unlink("splat.gp");
+ unlink("profile.gp");
+ unlink("reference.gp");
+ fprintf(stdout," Done!\n");
+ fflush(stdout);
+ }
+
+ else
+ fprintf(stderr,"\n*** ERROR: Error occurred invoking gnuplot!\n");
+}
+
+void GraphLongley(struct site source, struct site destination, char *name)
+{
+ /* This function invokes gnuplot to generate an appropriate
+ output file indicating the Longley-Rice model loss between
+ the source and destination locations. "filename" is the
+ name assigned to the output file generated by gnuplot.
+ The filename extension is used to set gnuplot's terminal
+ setting and output file type. If no extension is found,
+ .gif is assumed. */
+
+ int x, y, z, errnum, errflag=0;
+ char filename[255], term[15], ext[15], strmode[100], report_name[80];
+ double maxloss=-100000.0, minloss=100000.0, loss, haavt, angle;
+ FILE *fd=NULL, *fd2=NULL;
+
+ sprintf(report_name,"%s-to-%s.lro",source.name,destination.name);
+
+ for (x=0; report_name[x]!=0; x++)
+ if (report_name[x]==32 || report_name[x]==17 || report_name[x]==92 || report_name[x]==42 || report_name[x]==47)
+ report_name[x]='_';
+
+ fd2=fopen(report_name,"w");
+
+ fprintf(fd2,"\n\t--==[ SPLAT! v%s Longley-Rice Model Path Loss Report ]==--\n\n",splat_version);
+ fprintf(fd2,"Analysis of RF path conditions between %s and %s:\n",source.name, destination.name);
+ fprintf(fd2,"\n-------------------------------------------------------------------------\n\n");
+ fprintf(fd2,"Transmitter site: %s\n",source.name);
+ fprintf(fd2,"Site location: %.4f North / %.4f West",source.lat, source.lon);
+ fprintf(fd2, " (%s N / ", dec2dms(source.lat));
+ fprintf(fd2, "%s W)\n", dec2dms(source.lon));
+ fprintf(fd2,"Ground elevation: %.2f feet AMSL\n",GetElevation(source));
+ fprintf(fd2,"Antenna height: %.2f feet AGL / %.2f feet AMSL\n",source.alt, source.alt+GetElevation(source));
+
+ haavt=haat(source);
+
+ if (haavt>-4999.0)
+ fprintf(fd2,"Antenna height above average terrain: %.2f feet\n",haavt);
+
+ fprintf(fd2,"Distance to %s: %.2f miles.\n",destination.name,Distance(source,destination));
+ fprintf(fd2,"Azimuth to %s: %.2f degrees.\n",destination.name,Azimuth(source,destination));
+
+ angle=ElevationAngle(source,destination);
+
+ if (angle>=0.0)
+ fprintf(fd2,"Angle of elevation between %s and %s: %+.4f degrees.\n",source.name,destination.name,angle);
+
+ if (angle<0.0)
+ fprintf(fd2,"Angle of depression between %s and %s: %+.4f degrees.\n",source.name,destination.name,angle);
+
+ fprintf(fd2,"\n-------------------------------------------------------------------------\n\n");
+
+ /* Receiver */
+
+ fprintf(fd2,"Receiver site: %s\n",destination.name);
+ fprintf(fd2,"Site location: %.4f North / %.4f West",destination.lat, destination.lon);
+ fprintf(fd2, " (%s N / ", dec2dms(destination.lat));
+ fprintf(fd2, "%s W)\n", dec2dms(destination.lon));
+ fprintf(fd2,"Ground elevation: %.2f feet AMSL\n",GetElevation(destination));
+ fprintf(fd2,"Antenna height: %.2f feet AGL / %.2f feet AMSL\n",destination.alt, destination.alt+GetElevation(destination));
+
+ haavt=haat(destination);
+
+ if (haavt>-4999.0)
+ fprintf(fd2,"Antenna height above average terrain: %.2f feet\n",haavt);
+
+ fprintf(fd2,"Distance to %s: %.2f miles.\n",source.name,Distance(source,destination));
+ fprintf(fd2,"Azimuth to %s: %.2f degrees.\n",source.name,Azimuth(destination,source));
+
+ angle=ElevationAngle(destination,source);
+
+ if (angle>=0.0)
+ fprintf(fd2,"Angle of elevation between %s and %s: %+.4f degrees.\n",destination.name,source.name,angle);
+
+ if (angle<0.0)
+ fprintf(fd2,"Angle of depression between %s and %s: %+.4f degrees.\n",destination.name,source.name,angle);
+
+ fprintf(fd2,"\n-------------------------------------------------------------------------\n\n");
+
+ fprintf(fd2,"Longley-Rice path calculation parameters used in this analysis:\n\n");
+ fprintf(fd2,"Earth's Dielectric Constant: %.3lf\n",LR.eps_dielect);
+ fprintf(fd2,"Earth's Conductivity: %.3lf\n",LR.sgm_conductivity);
+ fprintf(fd2,"Atmospheric Bending Constant (N): %.3lf\n",LR.eno_ns_surfref);
+ fprintf(fd2,"Frequency: %.3lf (MHz)\n",LR.frq_mhz);
+ fprintf(fd2,"Radio Climate: %d (",LR.radio_climate);
+
+ switch (LR.radio_climate)
+ {
+ case 1:
+ fprintf(fd2,"Equatorial");
+ break;
+
+ case 2:
+ fprintf(fd2,"Continental Subtropical");
+ break;
+
+ case 3:
+ fprintf(fd2,"Maritime Subtropical");
+ break;
+
+ case 4:
+ fprintf(fd2,"Desert");
+ break;
+
+ case 5:
+ fprintf(fd2,"Continental Temperate");
+ break;
+
+ case 6:
+ fprintf(fd2,"Martitime Temperate, Over Land");
+ break;
+
+ case 7:
+ fprintf(fd2,"Maritime Temperate, Over Sea");
+ break;
+
+ default:
+ fprintf(fd2,"Unknown");
+ }
+
+ fprintf(fd2,")\nPolarization: %d (",LR.pol);
+
+ if (LR.pol==0)
+ fprintf(fd2,"Horizontal");
+
+ if (LR.pol==1)
+ fprintf(fd2,"Vertical");
+
+ fprintf(fd2,")\nFraction of Situations: %.1lf%c\n",LR.conf*100.0,37);
+ fprintf(fd2,"Fraction of Time: %.1lf%c\n",LR.rel*100.0,37);
+
+ fprintf(fd2,"\n-------------------------------------------------------------------------\n\n");
+
+ fprintf(fd2,"Analysis Results:\n\n");
+
+ ReadPath(source, destination); /* destination=RX, source=TX */
+
+ elev_l[1]=0.04*METERS_PER_MILE;
+
+ for (x=1; x<path.length; x++)
+ elev_l[x+1]=path.elevation[x]*METERS_PER_FOOT;
+
+ fprintf(fd2,"Distance (mi)\tLoss (dB)\tErrnum\tComment\n\n");
+
+ fd=fopen("profile.gp","w");
+
+ for (x=2; x<path.length; x++)
+ {
+ elev_l[0]=x-1;
+
+ point_to_point(elev_l, source.alt*METERS_PER_FOOT,
+ destination.alt*METERS_PER_FOOT,
+ LR.eps_dielect, LR.sgm_conductivity, LR.eno_ns_surfref,
+ LR.frq_mhz, LR.radio_climate, LR.pol, LR.conf, LR.rel,
+ loss, strmode, errnum);
+
+ /* Note: PASS BY REFERENCE ... loss and errnum are pass
+ by reference, only used in this file by this function */
+
+
+ fprintf(fd,"%f\t%f\n",path.distance[path.length-1]-path.distance[x],loss);
+ fprintf(fd2,"%7.2f\t\t%7.2f\t\t %d\t%s\n",path.distance[x],loss, errnum, strmode);
+ errflag|=errnum;
+
+ if (loss>maxloss)
+ maxloss=loss;
+
+ if (loss<minloss)
+ minloss=loss;
+ }
+
+ fclose(fd);
+
+ if (errflag)
+ {
+ fprintf(fd2,"\nNotes on \"errnum\"...\n\n");
+ fprintf(fd2," 0: No error. :-)\n");
+ fprintf(fd2," 1: Warning! Some parameters are nearly out of range.\n");
+ fprintf(fd2," Results should be used with caution.\n");
+ fprintf(fd2," 2: Note: Default parameters have been substituted for impossible ones.\n");
+ fprintf(fd2," 3: Warning! A combination of parameters is out of range.\n");
+ fprintf(fd2," Results are probably invalid.\n");
+ fprintf(fd2," Other: Warning! Some parameters are out of range.\n");
+ fprintf(fd2," Results are probably invalid.\n\nEnd of Report\n");
+ }
+
+ fprintf(stdout,"Longley-Rice Path Loss between %s and %s is %.2f db\n",source.name, destination.name, loss);
+ fprintf(stdout,"Path Loss Report written to: \"%s\"\n",report_name);
+ fflush(stdout);
+
+ fclose(fd2);
+
+ if (name[0]==0)
+ {
+ /* Default filename and output file type */
+
+ strncpy(filename,"loss\0",5);
+ strncpy(term,"gif\0",4);
+ strncpy(ext,"gif\0",4);
+ }
+
+ else
+ {
+ /* Grab extension and terminal type from "name" */
+
+ for (x=0; name[x]!='.' && name[x]!=0 && x<254; x++)
+ filename[x]=name[x];
+
+ if (name[x]=='.')
+ {
+ for (y=0, z=x, x++; name[x]!=0 && x<254 && y<14; x++, y++)
+ {
+ term[y]=tolower(name[x]);
+ ext[y]=term[y];
+ }
+
+ ext[y]=0;
+ term[y]=0;
+ filename[z]=0;
+ }
+
+ else
+ { /* No extension -- Default is gif */
+
+ filename[x]=0;
+ strncpy(term,"gif\0",4);
+ strncpy(ext,"gif\0",4);
+ }
+ }
+
+ /* Either .ps or .postscript may be used
+ as an extension for postscript output. */
+
+ if (strncmp(term,"postscript",10)==0)
+ strncpy(ext,"ps\0",3);
+
+ else if (strncmp(ext,"ps",2)==0)
+ strncpy(term,"postscript\0",11);
+
+ fprintf(stdout,"Writing \"%s.%s\"...",filename,ext);
+ fflush(stdout);
+
+ fd=fopen("splat.gp","w");
+
+ fprintf(fd,"set grid\n");
+ fprintf(fd,"set yrange [%2.3f to %2.3f]\n", minloss, maxloss);
+ fprintf(fd,"set term %s\n",term);
+ fprintf(fd,"set title \"SPLAT! Loss Profile\"\n");
+ fprintf(fd,"set xlabel \"Distance Between %s and %s (miles)\"\n",destination.name,source.name);
+ fprintf(fd,"set ylabel \"Longley-Rice Loss (dB)\"\n");
+ fprintf(fd,"set output \"%s.%s\"\n",filename,ext);
+ fprintf(fd,"plot \"profile.gp\" title \"Longley-Rice Loss\" with lines\n");
+
+ fclose(fd);
+
+ x=system("gnuplot splat.gp");
+
+ if (x!=-1)
+ {
+ unlink("splat.gp");
+ unlink("profile.gp");
+ unlink("reference.gp");
+
+ fprintf(stdout," Done!\n");
+ fflush(stdout);
+ }
+
+ else
+ fprintf(stderr,"\n*** ERROR: Error occurred invoking gnuplot!\n");
+}
+
+void ObstructionReport(struct site xmtr, struct site rcvr, char report)
+{
+ struct site result, result2, new_site;
+ double angle, haavt;
+ unsigned char block;
+ char report_name[80], string[255];
+ int x;
+ FILE *fd;
+
+ sprintf(report_name,"%s-to-%s.txt",xmtr.name,rcvr.name);
+
+ for (x=0; report_name[x]!=0; x++)
+ if (report_name[x]==32 || report_name[x]==17 || report_name[x]==92 || report_name[x]==42 || report_name[x]==47)
+ report_name[x]='_';
+
+ fd=fopen(report_name,"w");
+
+ fprintf(fd,"\n\t\t--==[ SPLAT! v%s Obstruction Report ]==--\n\n",splat_version);
+ fprintf(fd,"Analysis of line-of-sight path conditions between %s and %s:\n",xmtr.name, rcvr.name);
+ fprintf(fd,"\n-------------------------------------------------------------------------\n\n");
+ fprintf(fd,"Transmitter site: %s\n",xmtr.name);
+ fprintf(fd,"Site location: %.4f North / %.4f West",xmtr.lat, xmtr.lon);
+ fprintf(fd, " (%s N / ", dec2dms(xmtr.lat));
+ fprintf(fd, "%s W)\n", dec2dms(xmtr.lon));
+ fprintf(fd,"Ground elevation: %.2f feet AMSL\n",GetElevation(xmtr));
+ fprintf(fd,"Antenna height: %.2f feet AGL / %.2f feet AMSL\n",xmtr.alt, xmtr.alt+GetElevation(xmtr));
+
+ haavt=haat(xmtr);
+
+ if (haavt>-4999.0)
+ fprintf(fd,"Antenna height above average terrain: %.2f feet\n",haavt);
+
+ fprintf(fd,"Distance to %s: %.2f miles.\n",rcvr.name,Distance(xmtr,rcvr));
+ fprintf(fd,"Azimuth to %s: %.2f degrees.\n",rcvr.name,Azimuth(xmtr,rcvr));
+
+ angle=ElevationAngle(xmtr,rcvr);
+
+ if (angle>=0.0)
+ fprintf(fd,"Angle of elevation between %s and %s: %+.4f degrees.\n",xmtr.name,rcvr.name,angle);
+
+ if (angle<0.0)
+ fprintf(fd,"Angle of depression between %s and %s: %+.4f degrees.\n",xmtr.name,rcvr.name,angle);
+
+ fprintf(fd,"\n-------------------------------------------------------------------------\n\n");
+
+ /* Receiver */
+
+ fprintf(fd,"Receiver site: %s\n",rcvr.name);
+ fprintf(fd,"Site location: %.4f North / %.4f West",rcvr.lat, rcvr.lon);
+ fprintf(fd, " (%s N / ", dec2dms(rcvr.lat));
+ fprintf(fd, "%s W)\n", dec2dms(rcvr.lon));
+ fprintf(fd,"Ground elevation: %.2f feet AMSL\n",GetElevation(rcvr));
+ fprintf(fd,"Antenna height: %.2f feet AGL / %.2f feet AMSL\n",rcvr.alt, rcvr.alt+GetElevation(rcvr));
+
+ haavt=haat(rcvr);
+
+ if (haavt>-4999.0)
+ fprintf(fd,"Antenna height above average terrain: %.2f feet\n",haavt);
+
+ fprintf(fd,"Distance to %s: %.2f miles.\n",xmtr.name,Distance(xmtr,rcvr));
+ fprintf(fd,"Azimuth to %s: %.2f degrees.\n",xmtr.name,Azimuth(rcvr,xmtr));
+
+ angle=ElevationAngle(rcvr,xmtr);
+
+ if (angle>=0.0)
+ fprintf(fd,"Angle of elevation between %s and %s: %+.4f degrees.\n",rcvr.name,xmtr.name,angle);
+
+ if (angle<0.0)
+ fprintf(fd,"Angle of depression between %s and %s: %+.4f degrees.\n",rcvr.name,xmtr.name,angle);
+
+ fprintf(fd,"\n-------------------------------------------------------------------------\n\n");
+
+ if (report=='y')
+ {
+ /* Write an Obstruction Report */
+
+ new_site=rcvr;
+ result=los(xmtr,rcvr);
+ result2=result;
+ result2.alt-=1;
+ block=result.name[0];
+
+ if (block=='*')
+ fprintf(fd,"SPLAT! detected obstructions at:\n\n");
+
+ while (block=='*')
+ {
+ if (result.lat!=result2.lat || result.lon!=result2.lon || result.alt!=result2.alt)
+ fprintf(fd,"\t%.4f N, %.4f W, %5.2f miles, %6.2f feet AMSL.\n",result.lat, result.lon, Distance(rcvr,result), result.alt);
+
+ result2=result;
+ new_site.alt+=1.0;
+
+ /* Can you hear me now? :-) */
+
+ result=los(xmtr,new_site);
+ block=result.name[0];
+ }
+
+ if (new_site.alt!=rcvr.alt)
+ sprintf(string,"\nAntenna at %s must be raised to at least %.2f feet AGL\nto clear all obstructions detected by SPLAT!\n\n",rcvr.name, new_site.alt);
+ else
+ sprintf(string,"\nNo obstructions due to terrain were detected by SPLAT!\n\n");
+ }
+
+ fprintf(fd,"%s",string);
+
+ fclose(fd);
+
+ /* Display LOS status to terminal */
+
+ fprintf(stdout,"%sObstruction report written to: \"%s\"\n",string,report_name);
+ fflush(stdout);
+}
+
+void SiteReport(struct site xmtr)
+{
+ char report_name[80];
+ double terrain;
+ int x, azi;
+ FILE *fd;
+
+ sprintf(report_name,"%s-site_report.txt",xmtr.name);
+
+ for (x=0; report_name[x]!=0; x++)
+ if (report_name[x]==32 || report_name[x]==17 || report_name[x]==92 || report_name[x]==42 || report_name[x]==47)
+ report_name[x]='_';
+
+ fd=fopen(report_name,"w");
+
+ fprintf(fd,"\n\t--==[ SPLAT! v%s Site Analysis Report For: %s ]==--\n\n",splat_version,xmtr.name);
+
+ fprintf(fd,"---------------------------------------------------------------------------\n\n");
+ fprintf(fd,"Site location: %.4f North / %.4f West",xmtr.lat, xmtr.lon);
+ fprintf(fd, " (%s N / ",dec2dms(xmtr.lat));
+ fprintf(fd, "%s W)\n",dec2dms(xmtr.lon));
+ fprintf(fd,"Ground elevation: %.2f feet AMSL\n",GetElevation(xmtr));
+ fprintf(fd,"Antenna height: %.2f feet AGL / %.2f feet AMSL\n",xmtr.alt, xmtr.alt+GetElevation(xmtr));
+
+ terrain=haat(xmtr);
+
+ if (terrain>-4999.0)
+ {
+ fprintf(fd,"Antenna height above average terrain: %.2f feet\n\n",terrain);
+
+ /* Display the average terrain between 2 and 10 miles
+ from the transmitter site at azimuths of 0, 45, 90,
+ 135, 180, 225, 270, and 315 degrees. */
+
+ for (azi=0; azi<=315; azi+=45)
+ {
+ fprintf(fd,"Average terrain at %3d degrees azimuth: ",azi);
+ terrain=AverageTerrain(xmtr,(double)azi,2.0,10.0);
+
+ if (terrain>-4999.0)
+ fprintf(fd,"%.2f feet AMSL\n",terrain);
+ else
+ fprintf(fd,"No terrain\n");
+ }
+ }
+
+ fprintf(fd,"\n---------------------------------------------------------------------------\n\n");
+ fclose(fd);
+ fprintf(stdout,"\nSite analysis report written to: \"%s\"\n",report_name);
+}
+
+int main(char argc, char *argv[])
+{
+
+ int x, y, z=0;
+ unsigned char rxlat, rxlon, txlat, txlon, min_lat,
+ min_lon, max_lat, max_lon,
+ coverage=0, LRmap=0,
+ ext[20], terrain_plot=0,
+ elevation_plot=0, height_plot=0,
+ longley_plot=0, cities=0, bfs=0, txsites=0,
+ count, west_min, west_max, north_min, north_max,
+ report='y';
+
+ char mapfile[255], header[80], city_file[5][255],
+ elevation_file[255], height_file[255],
+ longley_file[255], terrain_file[255],
+ string[255], rxfile[255],
+ txfile[255], map=0, boundary_file[5][255];
+
+ double altitude=0.0, altitudeLR=0.0, tx_range=0.0,
+ rx_range=0.0, deg_range=0.0, deg_limit,
+ deg_range_lon, er_mult;
+ struct site tx_site[4], rx_site;
+ FILE *fd;
+
+ sprintf(header,"\n --==[ SPLAT! v%s Terrain Analysis Software (c) 1997-2004 KD2BD ]==--\n\n", splat_version);
+
+ if (argc==1)
+ {
+ fprintf(stdout, "%sAvailable Options...\n\n\t-t txsite(s).qth (max of 4)\n\t-r rxsite.qth\n",header);
+ fprintf(stdout,"\t-c plot coverage area(s) of TX(s) based on an RX antenna at X feet AGL\n");
+ fprintf(stdout,"\t-L plot path loss map of TX based on an RX antenna at X feet AGL\n");
+ fprintf(stdout,"\t-s filename(s) of city/site file(s) to import (max of 5)\n");
+ fprintf(stdout,"\t-b filename(s) of cartographic boundary file(s) to import (max of 5)\n");
+ fprintf(stdout,"\t-p filename of terrain profile graph to plot\n");
+ fprintf(stdout,"\t-e filename of terrain elevation graph to plot\n");
+ fprintf(stdout,"\t-h filename of terrain height graph to plot\n");
+ fprintf(stdout,"\t-l filename of Longley-Rice graph to plot\n");
+ fprintf(stdout,"\t-o filename of topographic map to generate (.ppm)\n");
+ fprintf(stdout,"\t-d sdf file directory path (overrides path in ~/.splat_path file)\n");
+ fprintf(stdout,"\t-n no analysis, brief report\n\t-N no analysis, no report\n");
+ fprintf(stdout,"\t-m earth radius multiplier\n");
+ fprintf(stdout,"\t-R modify default range for -c or -L (miles)\n\n");
+
+ fprintf(stdout,"Type 'man splat', or see the documentation for more details.\n\n");
+ fflush(stdout);
+ return 1;
+ }
+
+ y=argc-1;
+
+ rxfile[0]=0;
+ txfile[0]=0;
+ string[0]=0;
+ mapfile[0]=0;
+ elevation_file[0]=0;
+ terrain_file[0]=0;
+ sdf_path[0]=0;
+ path.length=0;
+ rx_site.lat=0.0;
+ rx_site.lon=0.0;
+ earthradius=EARTHRADIUS;
+
+ for (x=0; x<MAXSLOTS; x++)
+ {
+ dem[x].min_el=0;
+ dem[x].max_el=0;
+ dem[x].min_north=0;
+ dem[x].max_north=0;
+ dem[x].min_west=0;
+ dem[x].max_west=0;
+ }
+
+ /* Scan for command line arguments */
+
+ for (x=1; x<=y; x++)
+ {
+ if (strcmp(argv[x],"-R")==0)
+ {
+ z=x+1;
+
+ if (z<=y && argv[z][0] && argv[z][0]!='-')
+ {
+ sscanf(argv[z],"%lf",&max_range);
+
+ if (max_range<0.0)
+ max_range=0.0;
+
+ if (max_range>1000.0)
+ max_range=1000.0;
+ }
+ }
+
+ if (strcmp(argv[x],"-m")==0)
+ {
+ z=x+1;
+
+ if (z<=y && argv[z][0] && argv[z][0]!='-')
+ {
+ sscanf(argv[z],"%lf",&er_mult);
+
+ if (er_mult<0.1)
+ er_mult=1.0;
+
+ if (er_mult>1.0e6)
+ er_mult=1.0e6;
+
+ earthradius*=er_mult;
+ }
+ }
+
+ if (strcmp(argv[x],"-o")==0)
+ {
+ z=x+1;
+
+ if (z<=y && argv[z][0] && argv[z][0]!='-')
+ strncpy(mapfile,argv[z],253);
+ map=1;
+ }
+
+ if (strcmp(argv[x],"-c")==0)
+ {
+ z=x+1;
+
+ if (z<=y && argv[z][0] && argv[z][0]!='-')
+ {
+ sscanf(argv[z],"%lf",&altitude);
+ coverage=1;
+ }
+ }
+
+ if (strcmp(argv[x],"-p")==0)
+ {
+ z=x+1;
+
+ if (z<=y && argv[z][0] && argv[z][0]!='-')
+ {
+ strncpy(terrain_file,argv[z],253);
+ terrain_plot=1;
+ }
+ }
+
+ if (strcmp(argv[x],"-e")==0)
+ {
+ z=x+1;
+
+ if (z<=y && argv[z][0] && argv[z][0]!='-')
+ {
+ strncpy(elevation_file,argv[z],253);
+ elevation_plot=1;
+ }
+ }
+
+ if (strcmp(argv[x],"-h")==0)
+ {
+ z=x+1;
+
+ if (z<=y && argv[z][0] && argv[z][0]!='-')
+ {
+ strncpy(height_file,argv[z],253);
+ height_plot=1;
+ }
+ }
+
+ if (strcmp(argv[x],"-n")==0)
+ {
+ if (z<=y && argv[z][0] && argv[z][0]!='-')
+ {
+ report='n';
+ map=1;
+ }
+ }
+
+ if (strcmp(argv[x],"-N")==0)
+ {
+ if (z<=y && argv[z][0] && argv[z][0]!='-');
+ {
+ report='N';
+ map=1;
+ }
+ }
+
+ if (strcmp(argv[x],"-d")==0)
+ {
+ z=x+1;
+
+ if (z<=y && argv[z][0] && argv[z][0]!='-')
+ strncpy(sdf_path,argv[z],253);
+ }
+
+ if (strcmp(argv[x],"-t")==0)
+ {
+ /* Read Transmitter Location */
+
+ z=x+1;
+
+ while (z<=y && argv[z][0] && argv[z][0]!='-' && txsites<4)
+ {
+ strncpy(txfile,argv[z],253);
+ tx_site[txsites]=LoadQTH(txfile);
+ txsites++;
+ z++;
+ }
+ z--;
+ }
+
+ if (strcmp(argv[x],"-L")==0)
+ {
+ z=x+1;
+
+ if (z<=y && argv[z][0] && argv[z][0]!='-')
+ {
+ sscanf(argv[z],"%lf",&altitudeLR);
+
+ if (coverage)
+ fprintf(stdout,"c and L are exclusive options, ignoring L.\n");
+ else
+ {
+ LRmap=1;
+ ReadLRParm(txfile);
+ }
+ }
+ }
+
+ if (strcmp(argv[x],"-l")==0)
+ {
+ z=x+1;
+
+ if (z<=y && argv[z][0] && argv[z][0]!='-')
+ {
+ strncpy(longley_file,argv[z],253);
+ longley_plot=1;
+ /* Doing this twice is harmless */
+ ReadLRParm(txfile);
+ }
+ }
+
+ if (strcmp(argv[x],"-r")==0)
+ {
+ /* Read Receiver Location */
+
+ z=x+1;
+
+ if (z<=y && argv[z][0] && argv[z][0]!='-')
+ {
+ strncpy(rxfile,argv[z],253);
+ rx_site=LoadQTH(rxfile);
+ }
+ }
+
+ if (strcmp(argv[x],"-s")==0)
+ {
+ /* Read city file(s) */
+
+ z=x+1;
+
+ while (z<=y && argv[z][0] && argv[z][0]!='-' && cities<5)
+ {
+ strncpy(city_file[cities],argv[z],253);
+ cities++;
+ z++;
+ }
+ z--;
+ }
+
+ if (strcmp(argv[x],"-b")==0)
+ {
+ /* Read Boundary File(s) */
+
+ z=x+1;
+
+ while (z<=y && argv[z][0] && argv[z][0]!='-' && bfs<5)
+ {
+ strncpy(boundary_file[bfs],argv[z],253);
+ bfs++;
+ z++;
+ }
+ z--;
+ }
+ }
+
+ /* Perform some error checking on the arguments
+ and switches parsed from the command-line.
+ If an error is encountered, print a message
+ and exit gracefully. */
+
+ if (txsites==0)
+ {
+ fprintf(stderr,"\n%c*** ERROR: No transmitter site(s) specified!\n\n",7);
+ exit (-1);
+ }
+
+ for (x=0, y=0; x<txsites; x++)
+ {
+ if (tx_site[x].lat==0.0 && tx_site[x].lon==0.0)
+ {
+ fprintf(stderr,"\n*** ERROR: Transmitter site #%d not found!",x+1);
+ y++;
+ }
+ }
+
+ if (y)
+ {
+ fprintf(stderr,"%c\n\n",7);
+ exit (-1);
+ }
+
+ if ((coverage+LRmap)==0 && rx_site.lat==0.0 && rx_site.lon==0.0)
+ {
+ fprintf(stderr,"\n%c*** ERROR: No receiver site found or specified!\n\n",7);
+ exit (-1);
+ }
+
+ /* No errors were detected. Whew! :-) */
+
+ /* If no SDF path was specified on the command line (-d), check
+ for a path specified in the $HOME/.splat_path file. If the
+ file is not found, then sdf_path[] remains NULL, and the
+ current working directory is assumed to contain the SDF
+ files. */
+
+ if (sdf_path[0]==0)
+ {
+ sprintf(string,"%s/.splat_path",getenv("HOME"));
+ fd=fopen(string,"r");
+
+ if (fd!=NULL)
+ {
+ fgets(string,253,fd);
+
+ /* Remove <CR> and/or <LF> from string */
+
+ for (x=0; string[x]!=13 && string[x]!=10 && string[x]!=0 && x<253; x++);
+ string[x]=0;
+
+ strncpy(sdf_path,string,253);
+
+ fclose(fd);
+ }
+ }
+
+ /* Ensure a trailing '/' is present in sdf_path */
+
+ if (sdf_path[0])
+ {
+ x=strlen(sdf_path);
+
+ if (sdf_path[x-1]!='/' && x!=0)
+ {
+ sdf_path[x]='/';
+ sdf_path[x+1]=0;
+ }
+ }
+
+ fprintf(stdout,"%s",header);
+ fflush(stdout);
+
+ x=0;
+ y=0;
+
+ min_lat=0;
+ max_lat=0;
+ min_lon=0;
+ max_lon=0;
+
+ rxlat=(unsigned char)floor(rx_site.lat);
+ rxlon=(unsigned char)floor(rx_site.lon);
+
+ if (rxlat!=0)
+ {
+ if (min_lat==0)
+ min_lat=rxlat;
+
+ else if (rxlat<min_lat)
+ min_lat=rxlat;
+ }
+
+ if (rxlon!=0)
+ {
+ if (min_lon==0)
+ min_lon=rxlon;
+
+ else if (rxlon<min_lon)
+ min_lon=rxlon;
+ }
+
+ if (rxlat>max_lat)
+ max_lat=rxlat;
+
+ if (rxlon>max_lon)
+ max_lon=rxlon;
+
+ for (y=0, z=0; z<txsites; z++)
+ {
+ txlat=(unsigned char)floor(tx_site[z].lat);
+ txlon=(unsigned char)floor(tx_site[z].lon);
+
+ if (txlat!=0)
+ {
+ if (min_lat==0)
+ min_lat=txlat;
+
+ else if (txlat<min_lat)
+ min_lat=txlat;
+ }
+
+ if (txlon!=0)
+ {
+ if (min_lon==0)
+ min_lon=txlon;
+
+ else if (txlon<min_lon)
+ min_lon=txlon;
+ }
+
+ if (txlat>max_lat)
+ max_lat=txlat;
+
+ if (txlon>max_lon)
+ max_lon=txlon;
+ }
+
+ if (min_lat!=0 && min_lon!=0 && max_lat!=0 && max_lon!=0)
+ {
+ for (y=min_lon; y<=max_lon; y++)
+ for (x=min_lat; x<=max_lat; x++)
+ {
+ sprintf(string,"%u:%u:%u:%u",x, x+1, y, y+1);
+ LoadSDF(string);
+ }
+ }
+
+ if (coverage)
+ {
+ for (z=0; z<txsites; z++)
+ {
+ /* "Ball park" estimates used to load any additional
+ SDF files required to conduct this analysis. */
+
+ tx_range=sqrt(1.5*(tx_site[z].alt+GetElevation(tx_site[z])));
+ rx_range=sqrt(1.5*altitude);
+
+ /* deg_range determines the maximum
+ amount of topo data we read */
+
+ deg_range=(tx_range+rx_range)/69.0;
+
+ /* max_range sets the maximum size of the
+ analysis. A small, non-zero amount can
+ be used to shrink the size of the analysis
+ and limit the amount of topo data read by
+ SPLAT! A very large number will only increase
+ the width of the analysis, not the size of
+ the map. */
+
+ if (max_range==0.0)
+ max_range=tx_range+rx_range;
+
+ if (max_range<(tx_range+rx_range))
+ deg_range=max_range/69.0;
+
+ /* Prevent the demand for a really wide coverage
+ from allocating more slots than are available
+ in memory. */
+
+ switch (MAXSLOTS)
+ {
+ case 2: deg_limit=0.25;
+ break;
+
+ case 4: deg_limit=0.5;
+ break;
+
+ case 9: deg_limit=1.0;
+ break;
+
+ case 16: deg_limit=2.0;
+ break;
+
+ case 25: deg_limit=3.0;
+ }
+
+ if (tx_site[z].lat<70.0)
+ deg_range_lon=deg_range/cos(deg2rad*tx_site[z].lat);
+ else
+ deg_range_lon=deg_range/cos(deg2rad*70.0);
+
+ /* Correct for squares in degrees not being square in miles */
+
+ if (deg_range>deg_limit)
+ deg_range=deg_limit;
+
+ if (deg_range_lon>deg_limit)
+ deg_range_lon=deg_limit;
+
+ north_min=(unsigned char)floor(tx_site[z].lat-deg_range);
+ north_max=(unsigned char)floor(tx_site[z].lat+deg_range);
+ west_min=(unsigned char)floor(tx_site[z].lon-deg_range_lon);
+ west_max=(unsigned char)floor(tx_site[z].lon+deg_range_lon);
+
+ if (min_lat==0)
+ min_lat=north_min;
+
+ else if (north_min<min_lat)
+ min_lat=north_min;
+
+ if (min_lon==0)
+ min_lon=west_min;
+
+ else if (west_min<min_lon)
+ min_lon=west_min;
+
+ if (north_max>max_lat)
+ max_lat=north_max;
+
+ if (west_max>max_lon)
+ max_lon=west_max;
+ }
+
+ if (min_lat!=0 && min_lon!=0 && max_lat!=0 && max_lon!=0)
+ {
+ for (y=min_lon; y<=max_lon; y++)
+ for (x=min_lat; x<=max_lat; x++)
+ {
+ sprintf(string,"%u:%u:%u:%u",x, x+1, y, y+1);
+ LoadSDF(string);
+ }
+ }
+ }
+
+ if (LRmap)
+ {
+ /* "Ball park" estimates used to load any additional
+ SDF files required to conduct this analysis. */
+
+ tx_range=sqrt(1.5*(tx_site[0].alt+GetElevation(tx_site[0])));
+ rx_range=sqrt(1.5*altitudeLR);
+
+ /**
+ tx_range=sqrt(5.0*tx_site[0].alt);
+ rx_range=sqrt(5.0*altitudeLR);
+ **/
+
+ /* deg_range determines the maximum
+ amount of topo data we read */
+
+ deg_range=(tx_range+rx_range)/69.0;
+
+ /* max_range sets the maximum size of the
+ analysis. A small, non-zero amount can
+ be used to shrink the size of the analysis
+ and limit the amount of topo data read by
+ SPLAT! A very large number will only increase
+ the width of the analysis, not the size of
+ the map. */
+
+ if (max_range==0.0)
+ max_range=tx_range+rx_range;
+
+ if (max_range<(tx_range+rx_range))
+ deg_range=max_range/69.0;
+
+ /* Prevent the demand for a really wide coverage
+ from allocating more slots than are available
+ in memory. */
+
+ switch (MAXSLOTS)
+ {
+ case 2: deg_limit=0.25;
+ break;
+
+ case 4: deg_limit=0.5;
+ break;
+
+ case 9: deg_limit=1.0;
+ break;
+
+ case 16: deg_limit=2.0;
+ break;
+
+ case 25: deg_limit=3.0;
+ }
+
+ if (tx_site[0].lat<70.0)
+ deg_range_lon=deg_range/cos(deg2rad*tx_site[0].lat);
+ else
+ deg_range_lon=deg_range/cos(deg2rad*70.0);
+
+ /* Correct for squares in degrees not being square in miles */
+
+ if (deg_range>deg_limit)
+ deg_range=deg_limit;
+
+ if (deg_range_lon>deg_limit)
+ deg_range_lon=deg_limit;
+
+ north_min=(unsigned char)floor(tx_site[0].lat-deg_range);
+ north_max=(unsigned char)floor(tx_site[0].lat+deg_range);
+ west_min=(unsigned char)floor(tx_site[0].lon-deg_range_lon);
+ west_max=(unsigned char)floor(tx_site[0].lon+deg_range_lon);
+
+ if (min_lat==0)
+ min_lat=north_min;
+
+ else if (north_min<min_lat)
+ min_lat=north_min;
+
+ if (min_lon==0)
+ min_lon=west_min;
+
+ else if (west_min<min_lon)
+ min_lon=west_min;
+
+ if (north_max>max_lat)
+ max_lat=north_max;
+
+ if (west_max>max_lon)
+ max_lon=west_max;
+
+ if (min_lat!=0 && min_lon!=0 && max_lat!=0 && max_lon!=0)
+ {
+ for (y=min_lon; y<=max_lon; y++)
+ for (x=min_lat; x<=max_lat; x++)
+ {
+ sprintf(string,"%u:%u:%u:%u",x, x+1, y, y+1);
+ LoadSDF(string);
+ }
+ }
+ }
+
+ if (mapfile[0])
+ map=1;
+
+ if (coverage)
+ {
+ for (x=0; x<txsites; x++)
+ {
+ PlotCoverage(tx_site[x],altitude);
+ PlaceMarker(tx_site[x]);
+
+ if (report!='N')
+ SiteReport(tx_site[x]);
+ }
+
+ map=1;
+ }
+
+ else if (LRmap)
+ {
+ PlotLRMap(tx_site[0],altitudeLR);
+ PlaceMarker(tx_site[0]);
+
+ if (report!='N')
+ SiteReport(tx_site[0]);
+
+ map=1;
+ }
+
+ else
+ {
+ PlaceMarker(rx_site);
+
+ for (x=0; x<txsites; x++)
+ {
+ PlaceMarker(tx_site[x]);
+
+ if (report=='y')
+ {
+ switch (x)
+ {
+ case 0:
+ PlotPath(tx_site[x],rx_site,1);
+ break;
+
+ case 1:
+ PlotPath(tx_site[x],rx_site,8);
+ break;
+
+ case 2:
+ PlotPath(tx_site[x],rx_site,16);
+ break;
+
+ case 3:
+ PlotPath(tx_site[x],rx_site,32);
+ }
+ }
+
+ if (report!='N')
+ ObstructionReport(tx_site[x],rx_site,report);
+ }
+ }
+
+ if (map)
+ {
+ if (bfs)
+ {
+ for (x=0; x<bfs; x++)
+ LoadBoundaries(boundary_file[x]);
+ }
+
+ if (cities)
+ {
+ for (x=0; x<cities; x++)
+ LoadCities(city_file[x]);
+ }
+
+ if (!LRmap)
+ WritePPM(mapfile);
+ else
+ WritePPMLR(mapfile);
+ }
+
+ if (terrain_plot)
+ {
+ if (txsites>1)
+ {
+ for (x=0; terrain_file[x]!='.' && terrain_file[x]!=0 && x<80; x++);
+
+ if (terrain_file[x]=='.') /* extension */
+ {
+ ext[0]='.';
+ for (y=1, z=x, x++; terrain_file[x]!=0 && x<253 && y<14; x++, y++)
+ ext[y]=terrain_file[x];
+
+ ext[y]=0;
+ terrain_file[z]=0;
+ }
+
+ else
+ {
+ ext[0]=0; /* No extension */
+ terrain_file[x]=0;
+ }
+
+ for (count=0; count<txsites; count++)
+ {
+ sprintf(string,"%s-%c%s%c",terrain_file,'1'+count,ext,0);
+ GraphTerrain(tx_site[count],rx_site,string);
+ }
+ }
+
+ else
+ GraphTerrain(tx_site[0],rx_site,terrain_file);
+ }
+
+ if (elevation_plot)
+ {
+ if (txsites>1)
+ {
+ for (x=0; elevation_file[x]!='.' && elevation_file[x]!=0 && x<80; x++);
+
+ if (elevation_file[x]=='.') /* extension */
+ {
+ ext[0]='.';
+ for (y=1, z=x, x++; elevation_file[x]!=0 && x<253 && y<14; x++, y++)
+ ext[y]=elevation_file[x];
+
+ ext[y]=0;
+ elevation_file[z]=0;
+ }
+
+ else
+ {
+ ext[0]=0; /* No extension */
+ elevation_file[x]=0;
+ }
+
+ for (count=0; count<txsites; count++)
+ {
+ sprintf(string,"%s-%c%s%c",elevation_file,'1'+count,ext,0);
+ GraphElevation(tx_site[count],rx_site,string);
+ }
+ }
+
+ else
+ GraphElevation(tx_site[0],rx_site,elevation_file);
+ }
+
+ if (height_plot)
+ {
+ if (txsites>1)
+ {
+ for (x=0; height_file[x]!='.' && height_file[x]!=0 && x<80; x++);
+
+ if (height_file[x]=='.') /* extension */
+ {
+ ext[0]='.';
+ for (y=1, z=x, x++; height_file[x]!=0 && x<253 && y<14; x++, y++)
+ ext[y]=height_file[x];
+
+ ext[y]=0;
+ height_file[z]=0;
+ }
+
+ else
+ {
+ ext[0]=0; /* No extension */
+ height_file[x]=0;
+ }
+
+ for (count=0; count<txsites; count++)
+ {
+ sprintf(string,"%s-%c%s%c",height_file,'1'+count,ext,0);
+ GraphHeight(tx_site[count],rx_site,string);
+ }
+ }
+
+ else
+ GraphHeight(tx_site[0],rx_site,height_file);
+ }
+
+ if (longley_plot)
+ {
+ if (txsites>1)
+ {
+ for (x=0; longley_file[x]!='.' && longley_file[x]!=0 && x<80; x++);
+
+ if (longley_file[x]=='.') /* extension */
+ {
+ ext[0]='.';
+ for (y=1, z=x, x++; longley_file[x]!=0 && x<253 && y<14; x++, y++)
+ ext[y]=longley_file[x];
+
+ ext[y]=0;
+ longley_file[z]=0;
+ }
+
+ else
+ {
+ ext[0]=0; /* No extension */
+ longley_file[x]=0;
+ }
+
+ for (count=0; count<txsites; count++)
+ {
+ sprintf(string,"%s-%c%s%c",longley_file,'1'+count,ext,0);
+ GraphLongley(tx_site[count],rx_site,string);
+ }
+ }
+
+ else
+ GraphLongley(tx_site[0],rx_site,longley_file);
+ }
+
+ return 0;
+}