1 /****************************************************************************
2 * SPLAT: A Terrain Analysis Program *
3 * Copyright John A. Magliacane, KD2BD 1997-2004 *
4 * Last update: 24-Jan-2004 *
5 *****************************************************************************
7 * This program is free software; you can redistribute it and/or modify it *
8 * under the terms of the GNU General Public License as published by the *
9 * Free Software Foundation; either version 2 of the License or any later *
12 * This program is distributed in the hope that it will useful, but WITHOUT *
13 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or *
14 * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License *
17 *****************************************************************************
19 * Extensively modified by J. D. McDonald in Jan. 2004 to include *
20 * the Longley-Rice propagation model using C++ code from NTIA/ITS. *
22 * See: http://elbert.its.bldrdoc.gov/itm.html *
24 *****************************************************************************
25 g++ -Wall -O3 -s -lm -lbz2 -fomit-frame-pointer itm.cpp splat.cpp -o splat
26 *****************************************************************************/
36 #include "smallfont.h"
40 #define BZBUFFER 65536
43 #define ARRAYSIZE 4950
47 #define ARRAYSIZE 10870
51 #define ARRAYSIZE 19240
55 #define ARRAYSIZE 30025
58 char string[255], sdf_path[255], opened=0, *splat_version={"1.1.0"};
60 double TWOPI=6.283185307179586, HALFPI=1.570796326794896,
61 PI=3.141592653589793, deg2rad=1.74532925199e-02,
62 EARTHRADIUS=20902230.97, METERS_PER_MILE=1609.344,
63 METERS_PER_FOOT=0.3048, earthradius, max_range=0.0;
65 int min_north=0, max_north=0, min_west=0, max_west=0,
66 max_elevation=0, min_elevation=0, bzerror;
68 struct site { double lat;
74 struct { float lat[ARRAYSIZE];
76 float elevation[ARRAYSIZE];
77 float distance[ARRAYSIZE];
81 static struct { int min_north;
87 short data[1200][1200];
88 unsigned char mask[1200][1200];
93 double sgm_conductivity;
94 double eno_ns_surfref;
102 double elev_l[ARRAYSIZE+10];
104 void point_to_point(double elev[], double tht_m, double rht_m,
105 double eps_dielect, double sgm_conductivity, double eno_ns_surfref,
106 double frq_mhz, int radio_climate, int pol, double conf,
107 double rel, double &dbloss, char *strmode, int &errnum);
109 double arccos(double x, double y)
111 /* This function implements the arc cosine function,
112 returning a value between 0 and TWOPI. */
125 char *dec2dms(double decimal)
127 /* Converts decimal degrees to degrees, minutes, seconds,
128 (DMS) and returns the result as a character string. */
130 int degrees, minutes, seconds;
149 sprintf(string,"%d%c %d\' %d\"", degrees, 176, minutes, seconds);
153 int OrMask(double lat, double lon, int value)
155 /* Lines, text, markings, and coverage areas are stored in a
156 mask that is combined with topology data when topographic
157 maps are generated by SPLAT!. This function sets bits in
158 the mask based on the latitude and longitude of the area
161 int x, y, indx, minlat, minlon;
164 minlat=(int)floor(lat);
165 minlon=(int)floor(lon);
167 for (indx=0, found=0; indx<MAXSLOTS && found==0;)
168 if (minlat==dem[indx].min_north && minlon==dem[indx].min_west)
175 x=(int)(1199.0*(lat-floor(lat)));
176 y=(int)(1199.0*(lon-floor(lon)));
178 dem[indx].mask[x][y]|=value;
179 return (dem[indx].mask[x][y]);
185 int GetMask(double lat, double lon)
187 /* This function returns the mask bits based on the latitude
188 and longitude given. */
190 return (OrMask(lat,lon,0));
193 double GetElevation(struct site location)
195 /* This function returns the elevation (in feet) of any location
196 represented by the digital elevation model data in memory.
197 Function returns -5000.0 for locations not found in memory. */
200 int x, y, indx, minlat, minlon;
205 minlat=(int)floor(location.lat);
206 minlon=(int)floor(location.lon);
208 x=(int)(1199.0*(location.lat-floor(location.lat)));
209 y=(int)(1199.0*(location.lon-floor(location.lon)));
211 for (indx=0, found=0; indx<MAXSLOTS && found==0; indx++)
213 if (minlat==dem[indx].min_north && minlon==dem[indx].min_west)
215 elevation=3.28084*dem[indx].data[x][y];
223 double Distance(struct site site1, struct site site2)
225 /* This function returns the great circle distance
226 in miles between any two site locations. */
228 double lat1, lon1, lat2, lon2, distance;
230 lat1=site1.lat*deg2rad;
231 lon1=site1.lon*deg2rad;
232 lat2=site2.lat*deg2rad;
233 lon2=site2.lon*deg2rad;
235 distance=3959.0*acos(sin(lat1)*sin(lat2)+cos(lat1)*cos(lat2)*cos((lon1)-(lon2)));
240 double Azimuth(struct site source, struct site destination)
242 /* This function returns the azimuth (in degrees) to the
243 destination as seen from the location of the source. */
245 double dest_lat, dest_lon, src_lat, src_lon,
246 beta, azimuth, diff, num, den, fraction;
248 dest_lat=destination.lat*deg2rad;
249 dest_lon=destination.lon*deg2rad;
251 src_lat=source.lat*deg2rad;
252 src_lon=source.lon*deg2rad;
254 /* Calculate Surface Distance */
256 beta=acos(sin(src_lat)*sin(dest_lat)+cos(src_lat)*cos(dest_lat)*cos(src_lon-dest_lon));
258 /* Calculate Azimuth */
260 num=sin(dest_lat)-(sin(src_lat)*cos(beta));
261 den=cos(src_lat)*sin(beta);
264 /* Trap potential problems in acos() due to rounding */
272 /* Calculate azimuth */
274 azimuth=acos(fraction);
276 /* Reference it to True North */
278 diff=dest_lon-src_lon;
287 azimuth=TWOPI-azimuth;
289 return (azimuth/deg2rad);
292 double ElevationAngle(struct site local, struct site remote)
294 /* This function returns the angle of elevation (in degrees)
295 of the remote location as seen from the local site.
296 A positive result represents an angle of elevation (uptilt),
297 while a negative result represents an angle of depression
298 (downtilt), as referenced to a normal to the center of
301 register double a, b, dx;
303 a=GetElevation(remote)+remote.alt+earthradius;
304 b=GetElevation(local)+local.alt+earthradius;
306 dx=5280.0*Distance(local,remote);
308 /* Apply the Law of Cosines */
310 return ((180.0*(acos(((b*b)+(dx*dx)-(a*a))/(2.0*b*dx)))/PI)-90.0);
313 void ReadPath(struct site source, struct site destination)
315 /* This function generates a sequence of latitude and
316 longitude positions between a source location and
317 a destination along a great circle path, and stores
318 elevation and distance information for points along
319 that path in the "path" structure for later use. */
321 double azimuth, distance, lat1, lon1, beta,
322 den, num, lat2, lon2, total_distance;
324 struct site tempsite;
328 lat1=source.lat*deg2rad;
329 lon1=source.lon*deg2rad;
330 azimuth=Azimuth(source,destination)*deg2rad;
331 total_distance=Distance(source,destination);
333 for (distance=0; distance<=total_distance; distance+=0.04)
335 beta=distance/3959.0;
336 lat2=asin(sin(lat1)*cos(beta)+cos(azimuth)*sin(beta)*cos(lat1));
337 num=cos(beta)-(sin(lat1)*sin(lat2));
338 den=cos(lat1)*cos(lat2);
340 if (azimuth==0.0 && (beta>HALFPI-lat1))
343 else if (azimuth==HALFPI && (beta>HALFPI+lat1))
346 else if (fabs(num/den)>1.0)
351 if ((PI-azimuth)>=0.0)
352 lon2=lon1-arccos(num,den);
354 lon2=lon1+arccos(num,den);
366 x1=(int)(1199.0*(lat2-floor(lat2)));
367 y1=(int)(1199.0*(lon2-floor(lon2)));
373 path.elevation[c]=GetElevation(tempsite);
374 path.distance[c]=distance;
378 /* Make sure exact destination point is recorded at path.length-1 */
380 x1=(int)(1199.0*(destination.lat-floor(destination.lat)));
381 y1=(int)(1199.0*(destination.lon-floor(destination.lon)));
383 path.lat[c]=destination.lat;
384 path.lon[c]=destination.lon;
385 path.elevation[c]=GetElevation(destination);
386 path.distance[c]=total_distance;
392 double AverageTerrain(struct site source, double azimuthx, double start_distance, double end_distance)
394 /* This function returns the average terrain calculated in
395 the direction of "azimuth" (degrees) between "start_distance"
396 and "end_distance" (miles) from the source location. If
397 the terrain is all water (non-critical error), -5000.0 is
398 returned. If not enough SDF data has been loaded into
399 memory to complete the survey (critical error), then
400 -9999.0 is returned. */
402 int c, samples, endpoint;
403 double beta, lat1, lon1, lat2, lon2, num, den, azimuth, terrain=0.0;
404 struct site destination;
406 lat1=source.lat*deg2rad;
407 lon1=source.lon*deg2rad;
409 /* Generate a path of elevations between the source
410 location and the remote location provided. */
412 beta=end_distance/3959.0;
414 azimuth=deg2rad*azimuthx;
416 lat2=asin(sin(lat1)*cos(beta)+cos(azimuth)*sin(beta)*cos(lat1));
417 num=cos(beta)-(sin(lat1)*sin(lat2));
418 den=cos(lat1)*cos(lat2);
420 if (azimuth==0.0 && (beta>HALFPI-lat1))
423 else if (azimuth==HALFPI && (beta>HALFPI+lat1))
426 else if (fabs(num/den)>1.0)
431 if ((PI-azimuth)>=0.0)
432 lon2=lon1-arccos(num,den);
434 lon2=lon1+arccos(num,den);
446 destination.lat=lat2;
447 destination.lon=lon2;
449 /* If SDF data is missing for the endpoint of
450 the radial, then the average terrain cannot
451 be accurately calculated. Return -9999.0 */
453 if (GetElevation(destination)<-4999.0)
457 ReadPath(source,destination);
459 endpoint=path.length;
461 /* Shrink the length of the radial if the
462 outermost portion is not over U.S. land. */
464 for (c=endpoint-1; c>=0 && path.elevation[c]==0.0; c--);
468 for (c=0, samples=0; c<endpoint; c++)
470 if (path.distance[c]>=start_distance)
472 terrain+=path.elevation[c];
478 terrain=-5000.0; /* No land */
480 terrain=(terrain/(double)samples);
486 double haat(struct site antenna)
488 /* This function returns the antenna's Height Above Average
489 Terrain (HAAT) based on FCC Part 73.313(d). If a critical
490 error occurs, such as a lack of SDF data to complete the
491 survey, -5000.0 is returned. */
495 double terrain, avg_terrain, haat, sum=0.0;
497 /* Calculate the average terrain between 2 and 10 miles
498 from the antenna site at azimuths of 0, 45, 90, 135,
499 180, 225, 270, and 315 degrees. */
501 for (c=0, azi=0; azi<=315 && error==0; azi+=45)
503 terrain=AverageTerrain(antenna, (double)azi, 2.0, 10.0);
505 if (terrain<-9998.0) /* SDF data is missing */
508 if (terrain>-4999.0) /* It's land, not water */
510 sum+=terrain; /* Sum of averages */
519 avg_terrain=(sum/(double)c);
520 haat=(antenna.alt+GetElevation(antenna))-avg_terrain;
525 void PlaceMarker(struct site location)
527 /* This function places text and marker data in the mask array
528 for illustration on topographic maps generated by SPLAT!.
529 By default, SPLAT! centers text information BELOW the marker,
530 but may move it above, to the left, or to the right of the
531 marker depending on how much room is available on the map,
532 or depending on whether the area is already occupied by
533 another marker or label. If no room or clear space is
534 available on the map to place the marker and its associated
535 text, then the marker and text are not written to the map. */
538 char ok2print, occupied;
539 double x, y, lat, lon, textx=0.0, texty=0.0, xmin, xmax,
540 ymin, ymax, p1, p3, p6, p8, p12, p16, p24, label_length;
549 if (lat<xmax && lat>xmin && lon<ymax && lon>ymin)
561 /* Is Marker Position Clear Of Text Or Other Markers? */
563 for (x=lat-p3; (x<=xmax && x>=xmin && x<=lat+p3); x+=p1)
564 for (y=lon-p3; (y<=ymax && y>=ymin && y<=lon+p3); y+=p1)
565 occupied|=(GetMask(x,y)&2);
569 /* Determine Where Text Can Be Positioned */
571 /* label_length=length in pixels.
572 Each character is 8 pixels wide. */
574 label_length=p1*(double)(strlen(location.name)<<3);
576 if (((lon+label_length)<=ymax) && (lon-label_length)>=ymin)
578 /* Default: Centered Text */
580 texty=lon+label_length/2.0;
584 /* Position Text Below The Marker */
591 /* Is This Position Clear Of
592 Text Or Other Markers? */
594 for (a=0, occupied=0; a<16; a++)
596 for (b=0; b<(int)strlen(location.name); b++)
597 for (c=0; c<8; c++, y-=p1)
598 occupied|=(GetMask(x,y)&2);
612 /* Position Text Above The Marker */
619 /* Is This Position Clear Of
620 Text Or Other Markers? */
622 for (a=0, occupied=0; a<16; a++)
624 for (b=0; b<(int)strlen(location.name); b++)
625 for (c=0; c<8; c++, y-=p1)
626 occupied|=(GetMask(x,y)&2);
641 if ((lon-label_length)>=ymin)
643 /* Position Text To The
644 Right Of The Marker */
652 /* Is This Position Clear Of
653 Text Or Other Markers? */
655 for (a=0, occupied=0; a<16; a++)
657 for (b=0; b<(int)strlen(location.name); b++)
658 for (c=0; c<8; c++, y-=p1)
659 occupied|=(GetMask(x,y)&2);
673 /* Position Text To The
674 Left Of The Marker */
677 texty=lon+p8+(label_length);
682 /* Is This Position Clear Of
683 Text Or Other Markers? */
685 for (a=0, occupied=0; a<16; a++)
687 for (b=0; b<(int)strlen(location.name); b++)
688 for (c=0; c<8; c++, y-=p1)
689 occupied|=(GetMask(x,y)&2);
702 /* textx and texty contain the latitude and longitude
703 coordinates that describe the placement of the text
706 if (ok2print && textx!=0.0 && texty!=0.0)
713 for (a=0; a<16 && ok2print; a++)
715 for (b=0; b<(int)strlen(location.name); b++)
717 byte=fontdata[16*(location.name[b])+a];
719 for (c=128; c>0; c=c>>1, y-=p1)
727 /* Draw Square Marker Centered
728 On Location Specified */
730 for (x=lat-p3; (x<=xmax && x>=xmin && x<=lat+p3); x+=p1)
731 for (y=lon-p3; (y<=ymax && y>=ymin && y<=lon+p3); y+=p1)
738 double ReadBearing(char *input)
740 /* This function takes numeric input in the form of a character
741 string, and returns an equivalent bearing in degrees as a
742 decimal number (double). The input may either be expressed
743 in decimal format (40.139722) or degree, minute, second
744 format (40 08 23). This function also safely handles
745 extra spaces found either leading, trailing, or
746 embedded within the numbers expressed in the
747 input string. Decimal seconds are permitted. */
751 int a, b, length, degrees, minutes;
754 /* Copy "input" to "string", and ignore any extra
755 spaces that might be present in the process. */
758 length=strlen(input);
760 for (a=0, b=0; a<length && a<18; a++)
762 if ((input[a]!=32 && input[a]!='\n') || (input[a]==32 && input[a+1]!=32 && b!=0))
771 /* Count number of spaces in the clean string. */
773 length=strlen(string);
775 for (a=0, b=0; a<length; a++)
779 if (b==0) /* Decimal Format (40.139722) */
780 sscanf(string,"%lf",&bearing);
782 if (b==2) /* Degree, Minute, Second Format (40 08 23) */
784 sscanf(string,"%d %d %lf",°rees, &minutes, &seconds);
785 bearing=(double)degrees+((double)minutes/60)+(seconds/3600);
788 /* Anything else returns a 0.0 */
790 if (bearing>360.0 || bearing<0.0)
796 struct site LoadQTH(char *filename)
798 /* This function reads SPLAT! .qth (site location) files.
799 The latitude and longitude may be expressed either in
800 decimal degrees, or in degree, minute, second format.
801 Antenna height is assumed to be expressed in feet above
802 ground level (AGL), unless followed by the letter 'M',
803 or 'm', or by the word "meters" or "Meters", in which
804 case meters is assumed, and is handled accordingly. */
807 char string[50], qthfile[255];
808 struct site tempsite;
811 for (x=0; filename[x]!='.' && filename[x]!=0 && x<250; x++)
812 qthfile[x]=filename[x];
825 fd=fopen(qthfile,"r");
832 /* Strip <CR> and/or <LF> from end of site name */
834 for (x=0; string[x]!=13 && string[x]!=10 && string[x]!=0; tempsite.name[x]=string[x], x++);
840 tempsite.lat=ReadBearing(string);
844 tempsite.lon=ReadBearing(string);
850 /* Remove <CR> and/or <LF> from antenna height string */
851 for (x=0; string[x]!=13 && string[x]!=10 && string[x]!=0; x++);
855 /* Antenna height may either be in feet or meters.
856 If the letter 'M' or 'm' is discovered in
857 the string, then this is an indication that
858 the value given is expressed in meters, and
859 must be converted to feet before exiting. */
861 for (x=0; string[x]!='M' && string[x]!='m' && string[x]!=0 && x<48; x++);
862 if (string[x]=='M' || string[x]=='m')
865 sscanf(string,"%lf",&tempsite.alt);
866 tempsite.alt*=3.28084;
872 sscanf(string,"%lf",&tempsite.alt);
879 int LoadSDF_SDF(char *name)
881 /* This function reads uncompressed SPLAT Data Files (.sdf)
882 containing digital elevation model data into memory.
883 Elevation data, maximum and minimum elevations, and
884 quadrangle limits are stored in the first available
887 int x, y, data, indx, minlat, minlon, maxlat, maxlon;
888 char found, free_slot=0, line[20], sdf_file[255], path_plus_name[255];
891 for (x=0; name[x]!='.' && name[x]!=0 && x<250; x++)
896 /* Parse filename for minimum latitude and longitude values */
898 sscanf(sdf_file,"%d:%d:%d:%d",&minlat,&maxlat,&minlon,&maxlon);
906 /* Is it already in memory? */
908 for (indx=0, found=0; indx<MAXSLOTS && found==0; indx++)
910 if (minlat!=0 && minlon!=0)
912 if (minlat==dem[indx].min_north && minlon==dem[indx].min_west && maxlat==dem[indx].max_north && maxlon==dem[indx].max_west)
917 /* Is room available to load it? */
921 for (indx=0, free_slot=0; indx<MAXSLOTS && free_slot==0; indx++)
922 if (dem[indx].max_north==0 && dem[indx].max_west==0)
928 if (free_slot && found==0 && indx>=0 && indx<MAXSLOTS && minlat!=0 && minlon!=0)
930 strncpy(path_plus_name,sdf_path,255);
931 strncat(path_plus_name,sdf_file,255);
933 fd=fopen(path_plus_name,"rb");
937 fprintf(stdout,"Loading \"%s\" into slot %d...",path_plus_name,indx+1);
941 sscanf(line,"%d",&dem[indx].max_west);
944 sscanf(line,"%d",&dem[indx].min_north);
947 sscanf(line,"%d",&dem[indx].min_west);
950 sscanf(line,"%d",&dem[indx].max_north);
952 for (x=0; x<1200; x++)
953 for (y=0; y<1200; y++)
956 sscanf(line,"%d",&data);
957 dem[indx].data[x][y]=data;
959 if (data>dem[indx].max_el)
960 dem[indx].max_el=data;
962 if (dem[indx].min_el==0)
963 dem[indx].min_el=data;
966 if (data<dem[indx].min_el)
967 dem[indx].min_el=data;
973 if (min_elevation==0)
974 min_elevation=dem[indx].min_el;
978 if (dem[indx].min_el<min_elevation)
979 min_elevation=dem[indx].min_el;
982 if (dem[indx].max_el>max_elevation)
983 max_elevation=dem[indx].max_el;
985 if (dem[indx].max_north>max_north)
986 max_north=dem[indx].max_north;
988 if (dem[indx].max_west>max_west)
989 max_west=dem[indx].max_west;
992 min_north=dem[indx].min_north;
995 if (dem[indx].min_north<min_north)
996 min_north=dem[indx].min_north;
1000 min_west=dem[indx].min_west;
1003 if (dem[indx].min_west<min_west)
1004 min_west=dem[indx].min_west;
1007 fprintf(stdout," Done!\n");
1020 char *BZfgets(BZFILE *bzfd, unsigned length)
1022 /* This function returns at most one less than 'length' number
1023 of characters from a bz2 compressed file whose file descriptor
1024 is pointed to by *bzfd. In operation, a buffer is filled with
1025 uncompressed data (size = BZBUFFER), which is then parsed
1026 and doled out as NULL terminated character strings every time
1027 this function is invoked. A NULL string indicates an EOF
1028 or error condition. */
1030 static int x, y, nBuf;
1031 static char buffer[BZBUFFER+1], output[BZBUFFER+1];
1034 if (opened!=1 && bzerror==BZ_OK)
1036 /* First time through. Initialize everything! */
1047 if (x==nBuf && bzerror!=BZ_STREAM_END && bzerror==BZ_OK && opened)
1049 /* Uncompress data into a static buffer */
1051 nBuf=BZ2_bzRead(&bzerror, bzfd, buffer, BZBUFFER);
1056 /* Build a string from buffer contents */
1058 output[y]=buffer[x];
1060 if (output[y]=='\n' || output[y]==0 || y==(int)length-1)
1079 int LoadSDF_BZ(char *name)
1081 /* This function reads .bz2 compressed SPLAT Data Files containing
1082 digital elevation model data into memory. Elevation data,
1083 maximum and minimum elevations, and quadrangle limits are
1084 stored in the first available dem[] structure. */
1086 int x, y, data, indx, minlat, minlon, maxlat, maxlon;
1087 char found, free_slot=0, sdf_file[255], path_plus_name[255];
1091 for (x=0; name[x]!='.' && name[x]!=0 && x<247; x++)
1092 sdf_file[x]=name[x];
1096 /* Parse sdf_file name for minimum latitude and longitude values */
1098 sscanf(sdf_file,"%d:%d:%d:%d",&minlat,&maxlat,&minlon,&maxlon);
1110 /* Is it already in memory? */
1112 for (indx=0, found=0; indx<MAXSLOTS && found==0; indx++)
1114 if (minlat!=0 && minlon!=0)
1116 if (minlat==dem[indx].min_north && minlon==dem[indx].min_west && maxlat==dem[indx].max_north && maxlon==dem[indx].max_west)
1121 /* Is room available to load it? */
1125 for (indx=0, free_slot=0; indx<MAXSLOTS && free_slot==0; indx++)
1126 if (dem[indx].max_north==0 && dem[indx].max_west==0)
1132 if (free_slot && found==0 && indx>=0 && indx<MAXSLOTS && minlat!=0 && minlon!=0)
1134 strncpy(path_plus_name,sdf_path,255);
1135 strncat(path_plus_name,sdf_file,255);
1137 fd=fopen(path_plus_name,"rb");
1138 bzfd=BZ2_bzReadOpen(&bzerror,fd,0,0,NULL,0);
1140 if (fd!=NULL && bzerror==BZ_OK)
1142 fprintf(stdout,"Loading \"%s\" into slot %d...",path_plus_name,indx+1);
1145 sscanf(BZfgets(bzfd,255),"%d",&dem[indx].max_west);
1146 sscanf(BZfgets(bzfd,255),"%d",&dem[indx].min_north);
1147 sscanf(BZfgets(bzfd,255),"%d",&dem[indx].min_west);
1148 sscanf(BZfgets(bzfd,255),"%d",&dem[indx].max_north);
1150 for (x=0; x<1200; x++)
1151 for (y=0; y<1200; y++)
1153 sscanf(BZfgets(bzfd,20),"%d",&data);
1154 dem[indx].data[x][y]=data;
1156 if (data>dem[indx].max_el)
1157 dem[indx].max_el=data;
1159 if (dem[indx].min_el==0)
1160 dem[indx].min_el=data;
1163 if (data<dem[indx].min_el)
1164 dem[indx].min_el=data;
1170 BZ2_bzReadClose(&bzerror,bzfd);
1172 if (min_elevation==0)
1173 min_elevation=dem[indx].min_el;
1177 if (dem[indx].min_el<min_elevation)
1178 min_elevation=dem[indx].min_el;
1181 if (dem[indx].max_el>max_elevation)
1182 max_elevation=dem[indx].max_el;
1184 if (dem[indx].max_north>max_north)
1185 max_north=dem[indx].max_north;
1187 if (dem[indx].max_west>max_west)
1188 max_west=dem[indx].max_west;
1191 min_north=dem[indx].min_north;
1194 if (dem[indx].min_north<min_north)
1195 min_north=dem[indx].min_north;
1199 min_west=dem[indx].min_west;
1202 if (dem[indx].min_west<min_west)
1203 min_west=dem[indx].min_west;
1206 fprintf(stdout," Done!\n");
1217 char LoadSDF(char *name)
1219 /* This function loads the requested SDF file from the filesystem.
1220 It first tries to invoke the LoadSDF_SDF() function to load an
1221 uncompressed SDF file (since uncompressed files load slightly
1222 faster). If that attempt fails, then it tries to load a
1223 compressed SDF file by invoking the LoadSDF_BZ() function.
1224 If that fails, then we can assume that no elevation data
1225 exists for the region requested, and that the region
1226 requested must be entirely over water. */
1228 int x, y, indx, minlat, minlon, maxlat, maxlon;
1229 char found, free_slot=0;
1230 int return_value=-1;
1232 /* Try to load an uncompressed SDF first. */
1234 return_value=LoadSDF_SDF(name);
1236 /* If that fails, try loading a compressed SDF. */
1238 if (return_value==0 || return_value==-1)
1239 return_value=LoadSDF_BZ(name);
1241 /* If neither format can be found, then assume the area is water. */
1243 if (return_value==0 || return_value==-1)
1245 /* Parse SDF name for minimum latitude and longitude values */
1247 sscanf(name,"%d:%d:%d:%d",&minlat,&maxlat,&minlon,&maxlon);
1249 /* Is it already in memory? */
1251 for (indx=0, found=0; indx<MAXSLOTS && found==0; indx++)
1253 if (minlat!=0 && minlon!=0)
1255 if (minlat==dem[indx].min_north && minlon==dem[indx].min_west && maxlat==dem[indx].max_north && maxlon==dem[indx].max_west)
1260 /* Is room available to load it? */
1264 for (indx=0, free_slot=0; indx<MAXSLOTS && free_slot==0; indx++)
1265 if (dem[indx].max_north==0 && dem[indx].max_west==0)
1271 if (free_slot && found==0 && indx>=0 && indx<MAXSLOTS && minlat!=0 && minlon!=0)
1273 fprintf(stdout,"Region \"%s\" assumed as sea-level into slot %d...",name,indx+1);
1276 dem[indx].max_west=maxlon;
1277 dem[indx].min_north=minlat;
1278 dem[indx].min_west=minlon;
1279 dem[indx].max_north=maxlat;
1281 for (x=0; x<1200; x++)
1282 for (y=0; y<1200; y++)
1284 dem[indx].data[x][y]=0;
1286 if (dem[indx].min_el>0)
1290 if (dem[indx].min_el<min_elevation)
1291 min_elevation=dem[indx].min_el;
1293 if (dem[indx].max_el>max_elevation)
1294 max_elevation=dem[indx].max_el;
1296 if (dem[indx].max_north>max_north)
1297 max_north=dem[indx].max_north;
1299 if (dem[indx].max_west>max_west)
1300 max_west=dem[indx].max_west;
1303 min_north=dem[indx].min_north;
1306 if (dem[indx].min_north<min_north)
1307 min_north=dem[indx].min_north;
1311 min_west=dem[indx].min_west;
1314 if (dem[indx].min_west<min_west)
1315 min_west=dem[indx].min_west;
1318 fprintf(stdout," Done!\n");
1325 return return_value;
1328 void LoadCities(char *filename)
1330 /* This function reads SPLAT! city/site files, and plots
1331 the locations and names of the cities and site locations
1332 read on topographic maps generated by SPLAT! */
1335 char input[80], str[3][80];
1336 struct site city_site;
1339 fd=fopen(filename,"r");
1345 fprintf(stdout,"Reading \"%s\"... ",filename);
1348 while (fd!=NULL && feof(fd)==0)
1350 /* Parse line for name, latitude, and longitude */
1352 for (x=0, y=0, z=0; x<78 && input[x]!=0 && z<3; x++)
1354 if (input[x]!=',' && y<78)
1368 strncpy(city_site.name,str[0],49);
1369 city_site.lat=ReadBearing(str[1]);
1370 city_site.lon=ReadBearing(str[2]);
1373 PlaceMarker(city_site);
1379 fprintf(stdout,"Done!\n");
1383 fprintf(stderr,"*** ERROR: \"%s\": not found!\n",filename);
1386 void LoadBoundaries(char *filename)
1388 /* This function reads Cartographic Boundary Files available from
1389 the U.S. Census Bureau, and plots the data contained in those
1390 files on the PPM Map generated by SPLAT!. Such files contain
1391 the coordinates that describe the boundaries of cities,
1392 counties, and states. */
1395 double lat0, lon0, lat1, lon1;
1397 struct site source, destination;
1400 fd=fopen(filename,"r");
1404 fgets(string,78,fd);
1406 fprintf(stdout,"Reading \"%s\"... ",filename);
1411 fgets(string,78,fd);
1412 sscanf(string,"%lf %lf", &lon0, &lat0);
1413 fgets(string,78,fd);
1417 sscanf(string,"%lf %lf", &lon1, &lat1);
1424 destination.lat=lat1;
1425 destination.lon=lon1;
1427 ReadPath(source,destination);
1429 for (x=0; x<path.length; x++)
1430 OrMask(path.lat[x],path.lon[x],4);
1435 fgets(string,78,fd);
1437 } while (strncmp(string,"END",3)!=0 && feof(fd)==0);
1439 fgets(string,78,fd);
1441 } while (strncmp(string,"END",3)!=0 && feof(fd)==0);
1445 fprintf(stdout,"Done!\n");
1449 fprintf(stderr,"*** ERROR: \"%s\": not found!\n",filename);
1452 void ReadLRParm(char *txsite_filename)
1454 /* This function reads Longley-Rice parameter data for the
1455 transmitter site. The file name is the same as the txsite,
1456 except the filename extension is .lrp. If the needed file
1457 is not found, then the file "splat.lrp" is read from the
1458 current working directory. Failure to load this file will
1459 result in the default parameters hard coded into this
1460 being used and written to "splat.lrp". */
1463 char filename[255], lookup[256], string[80];
1465 FILE *fd=NULL, *outfile=NULL;
1467 /* Default parameters in case things go bad */
1469 LR.eps_dielect=15.0;
1470 LR.sgm_conductivity=0.005;
1471 LR.eno_ns_surfref=301.0;
1478 /* Modify txsite filename to one with a .lrp extension. */
1480 strncpy(filename,txsite_filename,255);
1482 for (x=0; filename[x]!='.' && filename[x]!=0 && filename[x]!='\n' && x<249; x++);
1490 /* Small lookup table to parse file, pass
1491 numeric data, and ignore comments. */
1493 for (x=0; x<=255; x++)
1496 /* Valid characters */
1498 for (x=48; x<=57; x++)
1504 fd=fopen(filename,"r");
1508 /* Load default "splat.lrp" file */
1510 strncpy(filename,"splat.lrp\0",10);
1511 fd=fopen(filename,"r");
1516 fgets(string,80,fd);
1518 for (x=0; lookup[(int)string[x]] && x<20; x++)
1519 string[x]=lookup[(int)string[x]];
1523 ok=sscanf(string,"%lf", &din);
1529 fgets(string,80,fd);
1531 for (x=0; lookup[(int)string[x]] && x<20; x++)
1532 string[x]=lookup[(int)string[x]];
1536 ok=sscanf(string,"%lf", &din);
1541 LR.sgm_conductivity=din;
1543 fgets(string,80,fd);
1545 for (x=0; lookup[(int)string[x]] && x<20; x++)
1546 string[x]=lookup[(int)string[x]];
1550 ok=sscanf(string,"%lf", &din);
1555 LR.eno_ns_surfref=din;
1557 fgets(string,80,fd);
1559 for (x=0; lookup[(int)string[x]] && x<20; x++)
1560 string[x]=lookup[(int)string[x]];
1564 ok=sscanf(string,"%lf", &din);
1571 fgets(string,80,fd);
1573 for (x=0; lookup[(int)string[x]] && x<20; x++)
1574 string[x]=lookup[(int)string[x]];
1578 ok=sscanf(string,"%d", &iin);
1583 LR.radio_climate=iin;
1585 fgets(string,80,fd);
1587 for (x=0; lookup[(int)string[x]] && x<20; x++)
1588 string[x]=lookup[(int)string[x]];
1592 ok=sscanf(string,"%d", &iin);
1599 fgets(string,80,fd);
1601 for (x=0; lookup[(int)string[x]] && x<20; x++)
1602 string[x]=lookup[(int)string[x]];
1606 ok=sscanf(string,"%lf", &din);
1613 fgets(string,80,fd);
1615 for (x=0; lookup[(int)string[x]] && x<20; x++)
1616 string[x]=lookup[(int)string[x]];
1620 ok=sscanf(string,"%lf", &din);
1631 /* Create a "splat.lrp" file since one
1632 could not be successfully loaded. */
1634 outfile=fopen("splat.lrp","w");
1636 fprintf(outfile,"%.3f\t; Earth Dielectric Constant (Relative permittivity)\n",LR.eps_dielect);
1638 fprintf(outfile,"%.3f\t; Earth Conductivity (Siemens per meter)\n", LR.sgm_conductivity);
1640 fprintf(outfile,"%.3f\t; Atmospheric Bending Constant (N-Units)\n",LR.eno_ns_surfref);
1642 fprintf(outfile,"%.3f\t; Frequency in MHz (20 MHz to 20 GHz)\n", LR.frq_mhz);
1644 fprintf(outfile,"%d\t; Radio Climate\n",LR.radio_climate);
1646 fprintf(outfile,"%d\t; Polarization (0 = Horizontal, 1 = Vertical)\n", LR.pol);
1648 fprintf(outfile,"%.2f\t; Fraction of situations\n",LR.conf);
1650 fprintf(outfile, "%.2f\t; Fraction of time\n",LR.rel);
1652 fprintf(outfile,"\nPlease consult SPLAT! documentation for the meaning and use of this data.\n");
1656 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);
1659 if (fd==NULL || ok==0)
1660 fprintf(stderr,"Longley-Rice default parameters have been assumed for this analysis.\n");
1663 struct site los(struct site source, struct site destination)
1665 /* This function determines whether a line-of-sight path
1666 unobstructed by terrain exists between source (transmitter)
1667 and destination (receiver) based on the geographical
1668 locations of the two sites, their respective antenna
1669 heights above ground, and the terrain between them.
1670 A site structure is returned upon completion. If the
1671 first character of site.name is ' ', then a clear path
1672 exists between source and destination. If the first
1673 character is '*', then an obstruction exists, and the
1674 site.lat and site.lon elements of the structure provide
1675 the geographical location of the obstruction. */
1679 struct site test, blockage;
1680 register double distance, tx_alt, rx_alt,
1681 cos_xmtr_angle, cos_test_angle, test_alt;
1683 ReadPath(source,destination);
1685 distance=5280.0*Distance(source,destination);
1686 tx_alt=earthradius+source.alt+GetElevation(source);
1687 rx_alt=earthradius+destination.alt+GetElevation(destination);
1689 /* Elevation angle of the xmtr (source) as seen by the rcvr */
1691 cos_xmtr_angle=((rx_alt*rx_alt)+(distance*distance)-(tx_alt*tx_alt))/(2.0*rx_alt*distance);
1693 /* Determine the elevation angle of each discrete location
1694 along the path between the receiver and transmitter.
1696 Since obstructions are more likely due to terrain effects
1697 closest to the receiver rather than farther away, we start
1698 looking for potential obstructions from the receiver's
1699 location, and work our way towards the transmitter.
1700 This loop is broken when the first obstruction is
1701 detected. If we can travel all the way to the transmitter
1702 without detecting an obstruction, then we have a clear
1703 unobstructed path between transmitter and receiver. */
1705 for (x=path.length-1, block=0; x>0 && block==0; x--)
1707 /* Build a structure for each test
1708 point along the path to be surveyed. */
1710 test.lat=path.lat[x];
1711 test.lon=path.lon[x];
1713 /* Measure the distance between the
1714 test point and the receiver locations */
1716 distance=5280.0*Distance(test,destination);
1717 test_alt=earthradius+path.elevation[x];
1719 /* Determine the cosine of the elevation of the test
1720 point as seen from the location of the receiver */
1722 cos_test_angle=((rx_alt*rx_alt)+(distance*distance)-(test_alt*test_alt))/(2.0*rx_alt*distance);
1724 /* If the elevation angle to the test point (as seen from
1725 the receiver) is greater than the elevation angle to the
1726 transmitter (as seen by the receiver), then we have a
1727 path obstructed by terrain. Note: Since we're comparing
1728 the cosines of these angles rather than the angles
1729 themselves (eliminating the call to acos() saves
1730 considerable time), the following "if" statement is
1731 reversed from what it would normally be if the angles
1734 if (cos_xmtr_angle>cos_test_angle)
1737 blockage.lat=path.lat[x];
1738 blockage.lon=path.lon[x];
1739 blockage.alt=path.elevation[x];
1740 blockage.name[0]='*';
1749 blockage.name[0]=' ';
1755 void PlotPath(struct site source, struct site destination, char mask_value)
1757 /* This function analyzes the path between the source and
1758 destination locations. It determines which points along
1759 the path have line-of-sight visibility to the source.
1760 Points along with path having line-of-sight visibility
1761 to the source at an AGL altitude equal to that of the
1762 destination location are stored by setting bit 1 in the
1763 mask[][] array, which are displayed in green when PPM
1764 maps are later generated by SPLAT!. */
1768 register double cos_xmtr_angle, cos_test_angle, test_alt;
1769 double distance, rx_alt, tx_alt;
1771 ReadPath(source,destination);
1773 for (y=0; y<path.length; y++)
1775 /* Test this point only if it hasn't been already
1776 tested and found to be free of obstructions. */
1778 if ((GetMask(path.lat[y],path.lon[y])&mask_value)==0)
1780 distance=5280.0*path.distance[y];
1781 tx_alt=earthradius+source.alt+path.elevation[0];
1782 rx_alt=earthradius+destination.alt+path.elevation[y];
1784 /* Calculate the cosine of the elevation of the
1785 transmitter as seen at the temp rx point. */
1787 cos_xmtr_angle=((rx_alt*rx_alt)+(distance*distance)-(tx_alt*tx_alt))/(2.0*rx_alt*distance);
1789 for (x=y, block=0; x>=0 && block==0; x--)
1791 distance=5280.0*(path.distance[y]-path.distance[x]);
1792 test_alt=earthradius+path.elevation[x];
1794 cos_test_angle=((rx_alt*rx_alt)+(distance*distance)-(test_alt*test_alt))/(2.0*rx_alt*distance);
1796 /* Compare these two angles to determine if
1797 a blockage exists. Since we're comparing
1798 the cosines of these angles rather than
1799 the angles themselves, the following "if"
1800 statement is reversed from what it would
1801 be if the actual angles were compared. */
1803 if (cos_xmtr_angle>cos_test_angle)
1808 OrMask(path.lat[y],path.lon[y],mask_value);
1813 void PlotLRPath(struct site source, struct site destination)
1815 /* This function plots the RF signal path loss
1816 between source and destination points based
1817 on the Longley-Rice propagation model. */
1823 ReadPath(source,destination);
1824 elev_l[1]=0.04*METERS_PER_MILE;
1826 for (x=0; x<path.length; x++)
1827 elev_l[x+2]=path.elevation[x]*METERS_PER_FOOT;
1829 for (y=0; y<path.length; y++)
1831 /* Test this point only if it
1832 has not already been tested. */
1834 if (GetMask(path.lat[y],path.lon[y])==0 && 0.04*y<=max_range)
1838 point_to_point(elev_l,source.alt*METERS_PER_FOOT,
1839 destination.alt*METERS_PER_FOOT,
1840 LR.eps_dielect, LR.sgm_conductivity, LR.eno_ns_surfref,
1841 LR.frq_mhz, LR.radio_climate, LR.pol, LR.conf, LR.rel,
1842 loss, strmode, errnum);
1844 /* Note: PASS BY REFERENCE ... loss and errnum are pass
1845 by reference, only used in this file by this function */
1857 OrMask(path.lat[y],path.lon[y],((unsigned char)(loss))<<3);
1860 else if (GetMask(path.lat[y],path.lon[y])==0 && 0.04*y>max_range)
1861 OrMask(path.lat[y],path.lon[y],1);
1865 void PlotCoverage(struct site source, double altitude)
1867 /* This function performs a 360 degree sweep around the
1868 transmitter site (source location), and plots the
1869 line-of-sight coverage of the transmitter on the SPLAT!
1870 generated topographic map based on a receiver located
1871 at the specified altitude (in feet AGL). Results
1872 are stored in memory, and written out in the form
1873 of a topographic map when the WritePPM() function
1874 is later invoked. */
1876 double lat, lon, one_pixel;
1877 static unsigned char mask_value;
1880 unsigned char symbol[4], x;
1882 /* Initialize mask_value */
1884 if (mask_value!=8 && mask_value!=16 && mask_value!=32)
1887 one_pixel=1.0/1200.0;
1896 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);
1899 /* 18.75=1200 pixels/degree divided by 64 loops
1900 per progress indicator symbol (.oOo) printed. */
1902 z=(int)(18.75*(max_west-min_west));
1904 for (lon=min_west, x=0; lon<=max_west; lon+=one_pixel)
1910 PlotPath(source,edge,mask_value);
1915 fprintf(stdout,"%c",symbol[x]);
1927 fprintf(stdout,"\n25%c to 50%c ",37,37);
1930 z=(int)(18.75*(max_north-min_north));
1932 for (lat=max_north, x=0; lat>=min_north; lat-=one_pixel)
1938 PlotPath(source,edge,mask_value);
1943 fprintf(stdout,"%c",symbol[x]);
1955 fprintf(stdout,"\n50%c to 75%c ",37,37);
1958 z=(int)(18.75*(max_west-min_west));
1960 for (lon=min_west, x=0; lon<=max_west; lon+=one_pixel)
1966 PlotPath(source,edge,mask_value);
1971 fprintf(stdout,"%c",symbol[x]);
1983 fprintf(stdout,"\n75%c to 100%c ",37,37);
1986 z=(int)(18.75*(max_north-min_north));
1988 for (lat=min_north, x=0; lat<=max_north; lat+=one_pixel)
1994 PlotPath(source,edge,mask_value);
1999 fprintf(stdout,"%c",symbol[x]);
2010 fprintf(stdout,"\nDone!\n");
2013 /* Assign next mask value */
2030 void PlotLRMap(struct site source, double altitude)
2032 /* This function performs a 360 degree sweep around the
2033 transmitter site (source location), and plots the
2034 Longley-Rice attenuation on the SPLAT! generated
2035 topographic map based on a receiver located at
2036 the specified altitude (in feet AGL). Results
2037 are stored in memory, and written out in the form
2038 of a topographic map when the WritePPMLR() function
2039 is later invoked. */
2043 double lat, lon, one_pixel;
2044 unsigned char symbol[4], x;
2046 one_pixel=1.0/1200.0;
2055 fprintf(stdout,"\nComputing Longley-Rice coverage of %s ", source.name);
2056 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);
2059 /* 18.75=1200 pixels/degree divided by 64 loops
2060 per progress indicator symbol (.oOo) printed. */
2062 z=(int)(18.75*(max_west-min_west));
2064 for (lon=min_west, x=0; lon<=max_west; lon+=one_pixel)
2070 PlotLRPath(source,edge);
2075 fprintf(stdout,"%c",symbol[x]);
2087 fprintf(stdout,"\n25%c to 50%c ",37,37);
2090 z=(int)(18.75*(max_north-min_north));
2092 for (lat=max_north, x=0; lat>=min_north; lat-=one_pixel)
2098 PlotLRPath(source,edge);
2103 fprintf(stdout,"%c",symbol[x]);
2115 fprintf(stdout,"\n50%c to 75%c ",37,37);
2118 z=(int)(18.75*(max_west-min_west));
2120 for (lon=min_west, x=0; lon<=max_west; lon+=one_pixel)
2126 PlotLRPath(source,edge);
2131 fprintf(stdout,"%c",symbol[x]);
2143 fprintf(stdout,"\n75%c to 100%c ",37,37);
2146 z=(int)(18.75*(max_north-min_north));
2148 for (lat=min_north, x=0; lat<=max_north; lat+=one_pixel)
2154 PlotLRPath(source,edge);
2159 fprintf(stdout,"%c",symbol[x]);
2170 fprintf(stdout,"\nDone!\n");
2174 void WritePPM(char *filename)
2176 /* This function generates a topographic map in Portable Pix Map
2177 (PPM) format based on logarithmically scaled topology data,
2178 as well as the content of flags held in the mask[][] array.
2179 The image created is rotated counter-clockwise 90 degrees
2180 from its representation in dem[][] so that north points
2181 up and east points right in the image generated. */
2183 int indx, x, x0, y0, minlat, minlon;
2184 unsigned width, height, output;
2185 unsigned char found, mask;
2187 double conversion, lat, lon, one_over_gamma, one_pixel;
2190 one_pixel=1.0/1200.0;
2191 one_over_gamma=1.0/GAMMA;
2192 conversion=255.0/pow((double)(max_elevation-min_elevation),one_over_gamma);
2194 width=1200*(max_west-min_west);
2195 height=1200*(max_north-min_north);
2198 strncpy(mapfile, "map.ppm\0",8);
2201 for (x=0; filename[x]!='.' && filename[x]!=0 && x<250; x++)
2202 mapfile[x]=filename[x];
2211 fd=fopen(mapfile,"wb");
2213 fprintf(fd,"P6\n%u %u\n255\n",width,height);
2215 fprintf(stdout,"\nWriting \"%s\" (%ux%u pixmap image)... ",mapfile,width,height);
2218 for (lat=(double)max_north; lat>=(double)min_north; lat-=one_pixel)
2220 for (lon=(double)max_west; lon>=(double)min_west; lon-=one_pixel)
2222 minlat=(int)floor(lat);
2223 minlon=(int)floor(lon);
2225 for (indx=0, found=0; indx<MAXSLOTS && found==0;)
2226 if (minlat==dem[indx].min_north && minlon==dem[indx].min_west)
2232 x0=(int)((1199.0*(lat-floor(lat)))+0.5);
2233 y0=(int)((1199.0*(lon-floor(lon)))+0.5);
2235 mask=dem[indx].mask[x0][y0];
2238 /* Text Labels: Red */
2239 fprintf(fd,"%c%c%c",255,0,0);
2242 /* County Boundaries: Light Cyan */
2243 fprintf(fd,"%c%c%c",128,128,255);
2245 else switch (mask&57)
2249 fprintf(fd,"%c%c%c",0,255,0);
2254 fprintf(fd,"%c%c%c",0,255,255);
2258 /* TX1 + TX2: Yellow */
2259 fprintf(fd,"%c%c%c",255,255,0);
2263 /* TX3: Medium Violet */
2264 fprintf(fd,"%c%c%c",147,112,219);
2268 /* TX1 + TX3: Pink */
2269 fprintf(fd,"%c%c%c",255,192,203);
2273 /* TX2 + TX3: Orange */
2274 fprintf(fd,"%c%c%c",255,165,0);
2278 /* TX1 + TX2 + TX3: Dark Green */
2279 fprintf(fd,"%c%c%c",0,100,0);
2284 fprintf(fd,"%c%c%c",255,130,71);
2288 /* TX1 + TX4: Green Yellow */
2289 fprintf(fd,"%c%c%c",173,255,47);
2293 /* TX2 + TX4: Dark Sea Green 1 */
2294 fprintf(fd,"%c%c%c",193,255,193);
2298 /* TX1 + TX2 + TX4: Blanched Almond */
2299 fprintf(fd,"%c%c%c",255,235,205);
2303 /* TX3 + TX4: Dark Turquoise */
2304 fprintf(fd,"%c%c%c",0,206,209);
2308 /* TX1 + TX3 + TX4: Medium Spring Green */
2309 fprintf(fd,"%c%c%c",0,250,154);
2313 /* TX2 + TX3 + TX4: Tan */
2314 fprintf(fd,"%c%c%c",210,180,140);
2318 /* TX1 + TX2 + TX3 + TX4: Gold2 */
2319 fprintf(fd,"%c%c%c",238,201,0);
2323 /* Water: Medium Blue */
2324 if (dem[indx].data[x0][y0]==0)
2325 fprintf(fd,"%c%c%c",0,0,170);
2328 /* Elevation: Greyscale */
2329 output=(unsigned)(0.5+pow((double)(dem[indx].data[x0][y0]-min_elevation),one_over_gamma)*conversion);
2330 fprintf(fd,"%c%c%c",output,output,output);
2338 fprintf(stdout,"Done!\n");
2342 void WritePPMLR(char *filename)
2344 /* This function generates a topographic map in Portable Pix Map
2345 (PPM) format based on the content of flags held in the mask[][]
2346 array (only). The image created is rotated counter-clockwise
2347 90 degrees from its representation in dem[][] so that north
2348 points up and east points right in the image generated. */
2350 int indx, x, t, t2, x0, y0, minlat, minlon;
2351 unsigned width, height, output;
2352 unsigned char found, mask;
2354 double conversion, lat, lon, one_over_gamma, one_pixel;
2357 one_pixel=1.0/1200.0;
2358 one_over_gamma=1.0/GAMMA;
2359 conversion=255.0/pow((double)(max_elevation-min_elevation),one_over_gamma);
2361 width=1200*(max_west-min_west);
2362 height=1200*(max_north-min_north);
2365 strncpy(mapfile, "map.ppm\0",8);
2368 for (x=0; filename[x]!='.' && filename[x]!=0 && x<250; x++)
2369 mapfile[x]=filename[x];
2378 fd=fopen(mapfile,"wb");
2380 fprintf(fd,"P6\n%u %u\n255\n",width,height+30);
2382 fprintf(stdout,"\nWriting \"%s\" (%ux%u pixmap image)... ",mapfile,width,height+30);
2385 for (lat=(double)max_north; lat>=(double)min_north; lat-=one_pixel)
2387 for (lon=(double)max_west; lon>=(double)min_west; lon-=one_pixel)
2389 minlat=(int)floor(lat);
2390 minlon=(int)floor(lon);
2392 for (indx=0, found=0; indx<MAXSLOTS && found==0;)
2393 if (minlat==dem[indx].min_north && minlon==dem[indx].min_west)
2399 x0=(int)((1199.0*(lat-floor(lat)))+0.5);
2400 y0=(int)((1199.0*(lon-floor(lon)))+0.5);
2402 mask=dem[indx].mask[x0][y0];
2406 /* Text Labels - Black or Red */
2408 fprintf(fd,"%c%c%c",0,0,0);
2410 fprintf(fd,"%c%c%c",255,0,0);
2414 /* County Boundaries: Black */
2415 fprintf(fd,"%c%c%c",0,0,0);
2417 else if (mask&1 && !((mask&248)==192))
2419 /* Outside Analysis Range */
2420 /* Display Greyscale / Sea Level */
2422 if (dem[indx].data[x0][y0]==0)
2423 fprintf(fd,"%c%c%c",0,0,170);
2426 output=(unsigned)(0.5+pow((double)(dem[indx].data[x0][y0]-min_elevation),one_over_gamma)*conversion);
2427 fprintf(fd,"%c%c%c",output,output,output);
2431 else switch ((mask&248)>>3)
2434 /* Inside range, but no coverage.
2435 Display Sea Level / Terrain */
2437 if (dem[indx].data[x0][y0]==0)
2438 fprintf(fd,"%c%c%c",0,0,170);
2441 /* Elevation: Greyscale */
2442 output=(unsigned)(0.5+pow((double)(dem[indx].data[x0][y0]-min_elevation),one_over_gamma)*conversion);
2443 fprintf(fd,"%c%c%c",output,output,output);
2450 fprintf(fd,"%c%c%c",0,255,0);
2455 fprintf(fd,"%c%c%c",255,192,203);
2460 fprintf(fd,"%c%c%c",0,255,255);
2465 fprintf(fd,"%c%c%c",255,255,0);
2470 fprintf(fd,"%c%c%c",161,131,224);
2475 fprintf(fd,"%c%c%c",255,165,0);
2480 fprintf(fd,"%c%c%c",193,255,193);
2485 fprintf(fd,"%c%c%c",255,108,108);
2489 /* TX1 + TX4: Green Yellow */
2490 fprintf(fd,"%c%c%c",173,255,47);
2494 /* Blanched Almond */
2495 fprintf(fd,"%c%c%c",255,235,184);
2499 /* Dark Turquoise */
2500 fprintf(fd,"%c%c%c",0,206,209);
2505 fprintf(fd,"%c%c%c",210,180,140);
2510 fprintf(fd,"%c%c%c",243,110,205);
2515 fprintf(fd,"%c%c%c",238,201,0);
2519 /* Medium Spring Green */
2520 fprintf(fd,"%c%c%c",0,250,154);
2524 /* Very light Blue */
2525 fprintf(fd,"%c%c%c",244,244,255);
2530 if (dem[indx].data[x0][y0]==0)
2531 fprintf(fd,"%c%c%c",0,0,170);
2534 /* Elevation: Greyscale */
2535 output=(unsigned)(0.5+pow((double)(dem[indx].data[x0][y0]-min_elevation),one_over_gamma)*conversion);
2536 fprintf(fd,"%c%c%c",output,output,output);
2545 for (y0=0; y0<30; y0++)
2547 for (indx=0; indx<16; indx++)
2549 for (x=0; x<x0; x++)
2554 if (y0>=10 && y0<=18)
2559 if (smallfont[t2/10][y0-10][x-11])
2564 if (smallfont[t2%10][y0-10][x-19])
2568 if (smallfont[0][y0-10][x-27])
2576 fprintf(fd,"%c%c%c",0,255,0);
2581 fprintf(fd,"%c%c%c",255,192,203);
2586 fprintf(fd,"%c%c%c",0,255,255);
2591 fprintf(fd,"%c%c%c",255,255,0);
2596 fprintf(fd,"%c%c%c",161,131,224);
2601 fprintf(fd,"%c%c%c",255,165,0);
2606 fprintf(fd,"%c%c%c",193,255,193);
2611 fprintf(fd,"%c%c%c",255,108,108);
2616 fprintf(fd,"%c%c%c",173,255,47);
2620 /* Blanched Almond */
2621 fprintf(fd,"%c%c%c",255,235,184);
2625 /* Dark Turquoise */
2626 fprintf(fd,"%c%c%c",0,206,209);
2631 fprintf(fd,"%c%c%c",210,180,140);
2636 fprintf(fd,"%c%c%c",243,110,205);
2641 fprintf(fd,"%c%c%c",238,201,0);
2645 /* Medium Spring Green */
2646 fprintf(fd,"%c%c%c",0,250,154);
2651 fprintf(fd,"%c%c%c",0,0,0);
2655 /* Very Light Blue */
2656 fprintf(fd,"%c%c%c",244,244,255);
2663 fprintf(stdout,"Done!\n");
2667 void GraphTerrain(struct site source, struct site destination, char *name)
2669 /* This function invokes gnuplot to generate an appropriate
2670 output file indicating the terrain profile between the source
2671 and destination locations. "filename" is the name assigned
2672 to the output file generated by gnuplot. The filename extension
2673 is used to set gnuplot's terminal setting and output file type.
2674 If no extension is found, .gif is assumed. */
2677 char filename[255], term[15], ext[15];
2680 ReadPath(destination,source);
2682 fd=fopen("profile.gp","wb");
2684 for (x=0; x<path.length; x++)
2685 fprintf(fd,"%f\t%f\n",path.distance[x],path.elevation[x]);
2691 /* Default filename and output file type */
2693 strncpy(filename,"profile\0",8);
2694 strncpy(term,"gif\0",4);
2695 strncpy(ext,"gif\0",4);
2700 /* Grab extension and terminal type from "name" */
2702 for (x=0; name[x]!='.' && name[x]!=0 && x<254; x++)
2703 filename[x]=name[x];
2707 for (y=0, z=x, x++; name[x]!=0 && x<254 && y<14; x++, y++)
2709 term[y]=tolower(name[x]);
2719 { /* No extension -- Default is gif */
2722 strncpy(term,"gif\0",4);
2723 strncpy(ext,"gif\0",4);
2727 /* Either .ps or .postscript may be used
2728 as an extension for postscript output. */
2730 if (strncmp(term,"postscript",10)==0)
2731 strncpy(ext,"ps\0",3);
2733 else if (strncmp(ext,"ps",2)==0)
2734 strncpy(term,"postscript\0",11);
2736 fprintf(stdout,"Writing \"%s.%s\"...",filename,ext);
2739 fd=fopen("splat.gp","w");
2740 fprintf(fd,"set grid\n");
2741 fprintf(fd,"set autoscale\n");
2742 fprintf(fd,"set term %s\n",term);
2743 fprintf(fd,"set title \"SPLAT! Terrain Profile\"\n");
2744 fprintf(fd,"set xlabel \"Distance Between %s and %s (miles)\"\n",destination.name,source.name);
2745 fprintf(fd,"set ylabel \"Ground Elevation Above Sea Level (feet)\"\n");
2746 fprintf(fd,"set output \"%s.%s\"\n",filename,ext);
2747 fprintf(fd,"plot \"profile.gp\" title \"\" with lines\n");
2750 x=system("gnuplot splat.gp");
2755 unlink("profile.gp");
2756 fprintf(stdout," Done!\n");
2761 fprintf(stderr,"\n*** ERROR: Error occurred invoking gnuplot!\n");
2764 void GraphElevation(struct site source, struct site destination, char *name)
2766 /* This function invokes gnuplot to generate an appropriate
2767 output file indicating the terrain profile between the source
2768 and destination locations. "filename" is the name assigned
2769 to the output file generated by gnuplot. The filename extension
2770 is used to set gnuplot's terminal setting and output file type.
2771 If no extension is found, .gif is assumed. */
2774 char filename[255], term[15], ext[15];
2775 double angle, refangle, maxangle=-90.0;
2777 FILE *fd=NULL, *fd2=NULL;
2779 ReadPath(destination,source); /* destination=RX, source=TX */
2780 refangle=ElevationAngle(destination,source);
2782 fd=fopen("profile.gp","wb");
2783 fd2=fopen("reference.gp","wb");
2785 for (x=1; x<path.length-1; x++)
2787 remote.lat=path.lat[x];
2788 remote.lon=path.lon[x];
2790 angle=ElevationAngle(destination,remote);
2791 fprintf(fd,"%f\t%f\n",path.distance[x],angle);
2792 fprintf(fd2,"%f\t%f\n",path.distance[x],refangle);
2798 fprintf(fd,"%f\t%f\n",path.distance[path.length-1],refangle);
2799 fprintf(fd2,"%f\t%f\n",path.distance[path.length-1],refangle);
2806 /* Default filename and output file type */
2808 strncpy(filename,"profile\0",8);
2809 strncpy(term,"gif\0",4);
2810 strncpy(ext,"gif\0",4);
2815 /* Grab extension and terminal type from "name" */
2817 for (x=0; name[x]!='.' && name[x]!=0 && x<254; x++)
2818 filename[x]=name[x];
2822 for (y=0, z=x, x++; name[x]!=0 && x<254 && y<14; x++, y++)
2824 term[y]=tolower(name[x]);
2834 { /* No extension -- Default is gif */
2837 strncpy(term,"gif\0",4);
2838 strncpy(ext,"gif\0",4);
2842 /* Either .ps or .postscript may be used
2843 as an extension for postscript output. */
2845 if (strncmp(term,"postscript",10)==0)
2846 strncpy(ext,"ps\0",3);
2848 else if (strncmp(ext,"ps",2)==0)
2849 strncpy(term,"postscript\0",11);
2851 fprintf(stdout,"Writing \"%s.%s\"...",filename,ext);
2854 fd=fopen("splat.gp","w");
2856 fprintf(fd,"set grid\n");
2857 fprintf(fd,"set yrange [%2.3f to %2.3f]\n", (-fabs(refangle)-0.25), maxangle+0.25);
2858 fprintf(fd,"set term %s\n",term);
2859 fprintf(fd,"set title \"SPLAT! Elevation Profile\"\n");
2860 fprintf(fd,"set xlabel \"Distance Between %s and %s (miles)\"\n",destination.name,source.name);
2861 fprintf(fd,"set ylabel \"Elevation Angle Along Path Between %s and %s (degrees)\"\n",destination.name,source.name);
2862 fprintf(fd,"set output \"%s.%s\"\n",filename,ext);
2863 fprintf(fd,"plot \"profile.gp\" title \"Real Earth Profile\" with lines, \"reference.gp\" title \"Line Of Sight Path\" with lines\n");
2867 x=system("gnuplot splat.gp");
2872 unlink("profile.gp");
2873 unlink("reference.gp");
2875 fprintf(stdout," Done!\n");
2880 fprintf(stderr,"\n*** ERROR: Error occurred invoking gnuplot!\n");
2883 void GraphHeight(struct site source, struct site destination, char *name)
2885 /* This function invokes gnuplot to generate an appropriate
2886 output file indicating the terrain profile between the source
2887 and destination locations. What is plotted is the height of land
2888 above or below a straight line between the receibe and transmit
2889 sites. "filename" is the name assigned to the output file
2890 generated by gnuplot. The filename extension is used to
2891 set gnuplot's terminal setting and output file type.
2892 If no extension is found, .gif is assumed. */
2895 char filename[255], term[15], ext[15];
2896 double a, b, c, height, refangle, cangle, maxheight=-100000.0,
2899 FILE *fd=NULL, *fd2=NULL;
2901 ReadPath(destination,source); /* destination=RX, source=TX */
2902 refangle=ElevationAngle(destination,source);
2903 b=GetElevation(destination)+destination.alt+earthradius;
2905 fd=fopen("profile.gp","wb");
2906 fd2=fopen("reference.gp","wb");
2908 for (x=1; x<path.length-1; x++)
2910 remote.lat=path.lat[x];
2911 remote.lon=path.lon[x];
2914 a=GetElevation(remote)+earthradius;
2916 cangle=5280.0*Distance(destination,remote)/earthradius;
2918 c=b*sin(refangle*deg2rad+HALFPI)/sin(HALFPI-refangle*deg2rad-cangle);
2922 fprintf(fd,"%f\t%f\n",path.distance[x],height);
2923 fprintf(fd2,"%f\t%f\n",path.distance[x],0.);
2925 if (height>maxheight)
2928 if (height<minheight)
2932 fprintf(fd,"%f\t%f\n",path.distance[path.length-1],0.0);
2933 fprintf(fd2,"%f\t%f\n",path.distance[path.length-1],0.0);
2940 /* Default filename and output file type */
2942 strncpy(filename,"height\0",8);
2943 strncpy(term,"gif\0",4);
2944 strncpy(ext,"gif\0",4);
2949 /* Grab extension and terminal type from "name" */
2951 for (x=0; name[x]!='.' && name[x]!=0 && x<254; x++)
2952 filename[x]=name[x];
2956 for (y=0, z=x, x++; name[x]!=0 && x<254 && y<14; x++, y++)
2958 term[y]=tolower(name[x]);
2968 { /* No extension -- Default is gif */
2971 strncpy(term,"gif\0",4);
2972 strncpy(ext,"gif\0",4);
2976 /* Either .ps or .postscript may be used
2977 as an extension for postscript output. */
2979 if (strncmp(term,"postscript",10)==0)
2980 strncpy(ext,"ps\0",3);
2982 else if (strncmp(ext,"ps",2)==0)
2983 strncpy(term,"postscript\0",11);
2985 fprintf(stdout,"Writing \"%s.%s\"...",filename,ext);
2988 fd=fopen("splat.gp","w");
2996 fprintf(fd,"set grid\n");
2997 fprintf(fd,"set yrange [%2.3f to %2.3f]\n", minheight, maxheight);
2998 fprintf(fd,"set term %s\n",term);
2999 fprintf(fd,"set title \"SPLAT! Height Profile\"\n");
3000 fprintf(fd,"set xlabel \"Distance Between %s and %s (miles)\"\n",destination.name,source.name);
3001 fprintf(fd,"set ylabel \"Ground Height Above Path Between %s and %s (feet)\"\n",destination.name,source.name);
3002 fprintf(fd,"set output \"%s.%s\"\n",filename,ext);
3003 fprintf(fd,"plot \"profile.gp\" title \"Real Earth Profile\" with lines, \"reference.gp\" title \"Line Of Sight Path\" with lines\n");
3007 x=system("gnuplot splat.gp");
3012 unlink("profile.gp");
3013 unlink("reference.gp");
3014 fprintf(stdout," Done!\n");
3019 fprintf(stderr,"\n*** ERROR: Error occurred invoking gnuplot!\n");
3022 void GraphLongley(struct site source, struct site destination, char *name)
3024 /* This function invokes gnuplot to generate an appropriate
3025 output file indicating the Longley-Rice model loss between
3026 the source and destination locations. "filename" is the
3027 name assigned to the output file generated by gnuplot.
3028 The filename extension is used to set gnuplot's terminal
3029 setting and output file type. If no extension is found,
3032 int x, y, z, errnum, errflag=0;
3033 char filename[255], term[15], ext[15], strmode[100], report_name[80];
3034 double maxloss=-100000.0, minloss=100000.0, loss, haavt, angle;
3035 FILE *fd=NULL, *fd2=NULL;
3037 sprintf(report_name,"%s-to-%s.lro",source.name,destination.name);
3039 for (x=0; report_name[x]!=0; x++)
3040 if (report_name[x]==32 || report_name[x]==17 || report_name[x]==92 || report_name[x]==42 || report_name[x]==47)
3043 fd2=fopen(report_name,"w");
3045 fprintf(fd2,"\n\t--==[ SPLAT! v%s Longley-Rice Model Path Loss Report ]==--\n\n",splat_version);
3046 fprintf(fd2,"Analysis of RF path conditions between %s and %s:\n",source.name, destination.name);
3047 fprintf(fd2,"\n-------------------------------------------------------------------------\n\n");
3048 fprintf(fd2,"Transmitter site: %s\n",source.name);
3049 fprintf(fd2,"Site location: %.4f North / %.4f West",source.lat, source.lon);
3050 fprintf(fd2, " (%s N / ", dec2dms(source.lat));
3051 fprintf(fd2, "%s W)\n", dec2dms(source.lon));
3052 fprintf(fd2,"Ground elevation: %.2f feet AMSL\n",GetElevation(source));
3053 fprintf(fd2,"Antenna height: %.2f feet AGL / %.2f feet AMSL\n",source.alt, source.alt+GetElevation(source));
3058 fprintf(fd2,"Antenna height above average terrain: %.2f feet\n",haavt);
3060 fprintf(fd2,"Distance to %s: %.2f miles.\n",destination.name,Distance(source,destination));
3061 fprintf(fd2,"Azimuth to %s: %.2f degrees.\n",destination.name,Azimuth(source,destination));
3063 angle=ElevationAngle(source,destination);
3066 fprintf(fd2,"Angle of elevation between %s and %s: %+.4f degrees.\n",source.name,destination.name,angle);
3069 fprintf(fd2,"Angle of depression between %s and %s: %+.4f degrees.\n",source.name,destination.name,angle);
3071 fprintf(fd2,"\n-------------------------------------------------------------------------\n\n");
3075 fprintf(fd2,"Receiver site: %s\n",destination.name);
3076 fprintf(fd2,"Site location: %.4f North / %.4f West",destination.lat, destination.lon);
3077 fprintf(fd2, " (%s N / ", dec2dms(destination.lat));
3078 fprintf(fd2, "%s W)\n", dec2dms(destination.lon));
3079 fprintf(fd2,"Ground elevation: %.2f feet AMSL\n",GetElevation(destination));
3080 fprintf(fd2,"Antenna height: %.2f feet AGL / %.2f feet AMSL\n",destination.alt, destination.alt+GetElevation(destination));
3082 haavt=haat(destination);
3085 fprintf(fd2,"Antenna height above average terrain: %.2f feet\n",haavt);
3087 fprintf(fd2,"Distance to %s: %.2f miles.\n",source.name,Distance(source,destination));
3088 fprintf(fd2,"Azimuth to %s: %.2f degrees.\n",source.name,Azimuth(destination,source));
3090 angle=ElevationAngle(destination,source);
3093 fprintf(fd2,"Angle of elevation between %s and %s: %+.4f degrees.\n",destination.name,source.name,angle);
3096 fprintf(fd2,"Angle of depression between %s and %s: %+.4f degrees.\n",destination.name,source.name,angle);
3098 fprintf(fd2,"\n-------------------------------------------------------------------------\n\n");
3100 fprintf(fd2,"Longley-Rice path calculation parameters used in this analysis:\n\n");
3101 fprintf(fd2,"Earth's Dielectric Constant: %.3lf\n",LR.eps_dielect);
3102 fprintf(fd2,"Earth's Conductivity: %.3lf\n",LR.sgm_conductivity);
3103 fprintf(fd2,"Atmospheric Bending Constant (N): %.3lf\n",LR.eno_ns_surfref);
3104 fprintf(fd2,"Frequency: %.3lf (MHz)\n",LR.frq_mhz);
3105 fprintf(fd2,"Radio Climate: %d (",LR.radio_climate);
3107 switch (LR.radio_climate)
3110 fprintf(fd2,"Equatorial");
3114 fprintf(fd2,"Continental Subtropical");
3118 fprintf(fd2,"Maritime Subtropical");
3122 fprintf(fd2,"Desert");
3126 fprintf(fd2,"Continental Temperate");
3130 fprintf(fd2,"Martitime Temperate, Over Land");
3134 fprintf(fd2,"Maritime Temperate, Over Sea");
3138 fprintf(fd2,"Unknown");
3141 fprintf(fd2,")\nPolarization: %d (",LR.pol);
3144 fprintf(fd2,"Horizontal");
3147 fprintf(fd2,"Vertical");
3149 fprintf(fd2,")\nFraction of Situations: %.1lf%c\n",LR.conf*100.0,37);
3150 fprintf(fd2,"Fraction of Time: %.1lf%c\n",LR.rel*100.0,37);
3152 fprintf(fd2,"\n-------------------------------------------------------------------------\n\n");
3154 fprintf(fd2,"Analysis Results:\n\n");
3156 ReadPath(source, destination); /* destination=RX, source=TX */
3158 elev_l[1]=0.04*METERS_PER_MILE;
3160 for (x=1; x<path.length; x++)
3161 elev_l[x+1]=path.elevation[x]*METERS_PER_FOOT;
3163 fprintf(fd2,"Distance (mi)\tLoss (dB)\tErrnum\tComment\n\n");
3165 fd=fopen("profile.gp","w");
3167 for (x=2; x<path.length; x++)
3171 point_to_point(elev_l, source.alt*METERS_PER_FOOT,
3172 destination.alt*METERS_PER_FOOT,
3173 LR.eps_dielect, LR.sgm_conductivity, LR.eno_ns_surfref,
3174 LR.frq_mhz, LR.radio_climate, LR.pol, LR.conf, LR.rel,
3175 loss, strmode, errnum);
3177 /* Note: PASS BY REFERENCE ... loss and errnum are pass
3178 by reference, only used in this file by this function */
3181 fprintf(fd,"%f\t%f\n",path.distance[path.length-1]-path.distance[x],loss);
3182 fprintf(fd2,"%7.2f\t\t%7.2f\t\t %d\t%s\n",path.distance[x],loss, errnum, strmode);
3196 fprintf(fd2,"\nNotes on \"errnum\"...\n\n");
3197 fprintf(fd2," 0: No error. :-)\n");
3198 fprintf(fd2," 1: Warning! Some parameters are nearly out of range.\n");
3199 fprintf(fd2," Results should be used with caution.\n");
3200 fprintf(fd2," 2: Note: Default parameters have been substituted for impossible ones.\n");
3201 fprintf(fd2," 3: Warning! A combination of parameters is out of range.\n");
3202 fprintf(fd2," Results are probably invalid.\n");
3203 fprintf(fd2," Other: Warning! Some parameters are out of range.\n");
3204 fprintf(fd2," Results are probably invalid.\n\nEnd of Report\n");
3207 fprintf(stdout,"Longley-Rice Path Loss between %s and %s is %.2f db\n",source.name, destination.name, loss);
3208 fprintf(stdout,"Path Loss Report written to: \"%s\"\n",report_name);
3215 /* Default filename and output file type */
3217 strncpy(filename,"loss\0",5);
3218 strncpy(term,"gif\0",4);
3219 strncpy(ext,"gif\0",4);
3224 /* Grab extension and terminal type from "name" */
3226 for (x=0; name[x]!='.' && name[x]!=0 && x<254; x++)
3227 filename[x]=name[x];
3231 for (y=0, z=x, x++; name[x]!=0 && x<254 && y<14; x++, y++)
3233 term[y]=tolower(name[x]);
3243 { /* No extension -- Default is gif */
3246 strncpy(term,"gif\0",4);
3247 strncpy(ext,"gif\0",4);
3251 /* Either .ps or .postscript may be used
3252 as an extension for postscript output. */
3254 if (strncmp(term,"postscript",10)==0)
3255 strncpy(ext,"ps\0",3);
3257 else if (strncmp(ext,"ps",2)==0)
3258 strncpy(term,"postscript\0",11);
3260 fprintf(stdout,"Writing \"%s.%s\"...",filename,ext);
3263 fd=fopen("splat.gp","w");
3265 fprintf(fd,"set grid\n");
3266 fprintf(fd,"set yrange [%2.3f to %2.3f]\n", minloss, maxloss);
3267 fprintf(fd,"set term %s\n",term);
3268 fprintf(fd,"set title \"SPLAT! Loss Profile\"\n");
3269 fprintf(fd,"set xlabel \"Distance Between %s and %s (miles)\"\n",destination.name,source.name);
3270 fprintf(fd,"set ylabel \"Longley-Rice Loss (dB)\"\n");
3271 fprintf(fd,"set output \"%s.%s\"\n",filename,ext);
3272 fprintf(fd,"plot \"profile.gp\" title \"Longley-Rice Loss\" with lines\n");
3276 x=system("gnuplot splat.gp");
3281 unlink("profile.gp");
3282 unlink("reference.gp");
3284 fprintf(stdout," Done!\n");
3289 fprintf(stderr,"\n*** ERROR: Error occurred invoking gnuplot!\n");
3292 void ObstructionReport(struct site xmtr, struct site rcvr, char report)
3294 struct site result, result2, new_site;
3295 double angle, haavt;
3296 unsigned char block;
3297 char report_name[80], string[255];
3301 sprintf(report_name,"%s-to-%s.txt",xmtr.name,rcvr.name);
3303 for (x=0; report_name[x]!=0; x++)
3304 if (report_name[x]==32 || report_name[x]==17 || report_name[x]==92 || report_name[x]==42 || report_name[x]==47)
3307 fd=fopen(report_name,"w");
3309 fprintf(fd,"\n\t\t--==[ SPLAT! v%s Obstruction Report ]==--\n\n",splat_version);
3310 fprintf(fd,"Analysis of line-of-sight path conditions between %s and %s:\n",xmtr.name, rcvr.name);
3311 fprintf(fd,"\n-------------------------------------------------------------------------\n\n");
3312 fprintf(fd,"Transmitter site: %s\n",xmtr.name);
3313 fprintf(fd,"Site location: %.4f North / %.4f West",xmtr.lat, xmtr.lon);
3314 fprintf(fd, " (%s N / ", dec2dms(xmtr.lat));
3315 fprintf(fd, "%s W)\n", dec2dms(xmtr.lon));
3316 fprintf(fd,"Ground elevation: %.2f feet AMSL\n",GetElevation(xmtr));
3317 fprintf(fd,"Antenna height: %.2f feet AGL / %.2f feet AMSL\n",xmtr.alt, xmtr.alt+GetElevation(xmtr));
3322 fprintf(fd,"Antenna height above average terrain: %.2f feet\n",haavt);
3324 fprintf(fd,"Distance to %s: %.2f miles.\n",rcvr.name,Distance(xmtr,rcvr));
3325 fprintf(fd,"Azimuth to %s: %.2f degrees.\n",rcvr.name,Azimuth(xmtr,rcvr));
3327 angle=ElevationAngle(xmtr,rcvr);
3330 fprintf(fd,"Angle of elevation between %s and %s: %+.4f degrees.\n",xmtr.name,rcvr.name,angle);
3333 fprintf(fd,"Angle of depression between %s and %s: %+.4f degrees.\n",xmtr.name,rcvr.name,angle);
3335 fprintf(fd,"\n-------------------------------------------------------------------------\n\n");
3339 fprintf(fd,"Receiver site: %s\n",rcvr.name);
3340 fprintf(fd,"Site location: %.4f North / %.4f West",rcvr.lat, rcvr.lon);
3341 fprintf(fd, " (%s N / ", dec2dms(rcvr.lat));
3342 fprintf(fd, "%s W)\n", dec2dms(rcvr.lon));
3343 fprintf(fd,"Ground elevation: %.2f feet AMSL\n",GetElevation(rcvr));
3344 fprintf(fd,"Antenna height: %.2f feet AGL / %.2f feet AMSL\n",rcvr.alt, rcvr.alt+GetElevation(rcvr));
3349 fprintf(fd,"Antenna height above average terrain: %.2f feet\n",haavt);
3351 fprintf(fd,"Distance to %s: %.2f miles.\n",xmtr.name,Distance(xmtr,rcvr));
3352 fprintf(fd,"Azimuth to %s: %.2f degrees.\n",xmtr.name,Azimuth(rcvr,xmtr));
3354 angle=ElevationAngle(rcvr,xmtr);
3357 fprintf(fd,"Angle of elevation between %s and %s: %+.4f degrees.\n",rcvr.name,xmtr.name,angle);
3360 fprintf(fd,"Angle of depression between %s and %s: %+.4f degrees.\n",rcvr.name,xmtr.name,angle);
3362 fprintf(fd,"\n-------------------------------------------------------------------------\n\n");
3366 /* Write an Obstruction Report */
3369 result=los(xmtr,rcvr);
3372 block=result.name[0];
3375 fprintf(fd,"SPLAT! detected obstructions at:\n\n");
3379 if (result.lat!=result2.lat || result.lon!=result2.lon || result.alt!=result2.alt)
3380 fprintf(fd,"\t%.4f N, %.4f W, %5.2f miles, %6.2f feet AMSL.\n",result.lat, result.lon, Distance(rcvr,result), result.alt);
3385 /* Can you hear me now? :-) */
3387 result=los(xmtr,new_site);
3388 block=result.name[0];
3391 if (new_site.alt!=rcvr.alt)
3392 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);
3394 sprintf(string,"\nNo obstructions due to terrain were detected by SPLAT!\n\n");
3397 fprintf(fd,"%s",string);
3401 /* Display LOS status to terminal */
3403 fprintf(stdout,"%sObstruction report written to: \"%s\"\n",string,report_name);
3407 void SiteReport(struct site xmtr)
3409 char report_name[80];
3414 sprintf(report_name,"%s-site_report.txt",xmtr.name);
3416 for (x=0; report_name[x]!=0; x++)
3417 if (report_name[x]==32 || report_name[x]==17 || report_name[x]==92 || report_name[x]==42 || report_name[x]==47)
3420 fd=fopen(report_name,"w");
3422 fprintf(fd,"\n\t--==[ SPLAT! v%s Site Analysis Report For: %s ]==--\n\n",splat_version,xmtr.name);
3424 fprintf(fd,"---------------------------------------------------------------------------\n\n");
3425 fprintf(fd,"Site location: %.4f North / %.4f West",xmtr.lat, xmtr.lon);
3426 fprintf(fd, " (%s N / ",dec2dms(xmtr.lat));
3427 fprintf(fd, "%s W)\n",dec2dms(xmtr.lon));
3428 fprintf(fd,"Ground elevation: %.2f feet AMSL\n",GetElevation(xmtr));
3429 fprintf(fd,"Antenna height: %.2f feet AGL / %.2f feet AMSL\n",xmtr.alt, xmtr.alt+GetElevation(xmtr));
3433 if (terrain>-4999.0)
3435 fprintf(fd,"Antenna height above average terrain: %.2f feet\n\n",terrain);
3437 /* Display the average terrain between 2 and 10 miles
3438 from the transmitter site at azimuths of 0, 45, 90,
3439 135, 180, 225, 270, and 315 degrees. */
3441 for (azi=0; azi<=315; azi+=45)
3443 fprintf(fd,"Average terrain at %3d degrees azimuth: ",azi);
3444 terrain=AverageTerrain(xmtr,(double)azi,2.0,10.0);
3446 if (terrain>-4999.0)
3447 fprintf(fd,"%.2f feet AMSL\n",terrain);
3449 fprintf(fd,"No terrain\n");
3453 fprintf(fd,"\n---------------------------------------------------------------------------\n\n");
3455 fprintf(stdout,"\nSite analysis report written to: \"%s\"\n",report_name);
3458 int main(char argc, char *argv[])
3462 unsigned char rxlat, rxlon, txlat, txlon, min_lat,
3463 min_lon, max_lat, max_lon,
3464 coverage=0, LRmap=0,
3465 ext[20], terrain_plot=0,
3466 elevation_plot=0, height_plot=0,
3467 longley_plot=0, cities=0, bfs=0, txsites=0,
3468 count, west_min, west_max, north_min, north_max,
3471 char mapfile[255], header[80], city_file[5][255],
3472 elevation_file[255], height_file[255],
3473 longley_file[255], terrain_file[255],
3474 string[255], rxfile[255],
3475 txfile[255], map=0, boundary_file[5][255];
3477 double altitude=0.0, altitudeLR=0.0, tx_range=0.0,
3478 rx_range=0.0, deg_range=0.0, deg_limit,
3479 deg_range_lon, er_mult;
3480 struct site tx_site[4], rx_site;
3483 sprintf(header,"\n --==[ SPLAT! v%s Terrain Analysis Software (c) 1997-2004 KD2BD ]==--\n\n", splat_version);
3487 fprintf(stdout, "%sAvailable Options...\n\n\t-t txsite(s).qth (max of 4)\n\t-r rxsite.qth\n",header);
3488 fprintf(stdout,"\t-c plot coverage area(s) of TX(s) based on an RX antenna at X feet AGL\n");
3489 fprintf(stdout,"\t-L plot path loss map of TX based on an RX antenna at X feet AGL\n");
3490 fprintf(stdout,"\t-s filename(s) of city/site file(s) to import (max of 5)\n");
3491 fprintf(stdout,"\t-b filename(s) of cartographic boundary file(s) to import (max of 5)\n");
3492 fprintf(stdout,"\t-p filename of terrain profile graph to plot\n");
3493 fprintf(stdout,"\t-e filename of terrain elevation graph to plot\n");
3494 fprintf(stdout,"\t-h filename of terrain height graph to plot\n");
3495 fprintf(stdout,"\t-l filename of Longley-Rice graph to plot\n");
3496 fprintf(stdout,"\t-o filename of topographic map to generate (.ppm)\n");
3497 fprintf(stdout,"\t-d sdf file directory path (overrides path in ~/.splat_path file)\n");
3498 fprintf(stdout,"\t-n no analysis, brief report\n\t-N no analysis, no report\n");
3499 fprintf(stdout,"\t-m earth radius multiplier\n");
3500 fprintf(stdout,"\t-R modify default range for -c or -L (miles)\n\n");
3502 fprintf(stdout,"Type 'man splat', or see the documentation for more details.\n\n");
3513 elevation_file[0]=0;
3519 earthradius=EARTHRADIUS;
3521 for (x=0; x<MAXSLOTS; x++)
3531 /* Scan for command line arguments */
3533 for (x=1; x<=y; x++)
3535 if (strcmp(argv[x],"-R")==0)
3539 if (z<=y && argv[z][0] && argv[z][0]!='-')
3541 sscanf(argv[z],"%lf",&max_range);
3546 if (max_range>1000.0)
3551 if (strcmp(argv[x],"-m")==0)
3555 if (z<=y && argv[z][0] && argv[z][0]!='-')
3557 sscanf(argv[z],"%lf",&er_mult);
3565 earthradius*=er_mult;
3569 if (strcmp(argv[x],"-o")==0)
3573 if (z<=y && argv[z][0] && argv[z][0]!='-')
3574 strncpy(mapfile,argv[z],253);
3578 if (strcmp(argv[x],"-c")==0)
3582 if (z<=y && argv[z][0] && argv[z][0]!='-')
3584 sscanf(argv[z],"%lf",&altitude);
3589 if (strcmp(argv[x],"-p")==0)
3593 if (z<=y && argv[z][0] && argv[z][0]!='-')
3595 strncpy(terrain_file,argv[z],253);
3600 if (strcmp(argv[x],"-e")==0)
3604 if (z<=y && argv[z][0] && argv[z][0]!='-')
3606 strncpy(elevation_file,argv[z],253);
3611 if (strcmp(argv[x],"-h")==0)
3615 if (z<=y && argv[z][0] && argv[z][0]!='-')
3617 strncpy(height_file,argv[z],253);
3622 if (strcmp(argv[x],"-n")==0)
3624 if (z<=y && argv[z][0] && argv[z][0]!='-')
3631 if (strcmp(argv[x],"-N")==0)
3633 if (z<=y && argv[z][0] && argv[z][0]!='-');
3640 if (strcmp(argv[x],"-d")==0)
3644 if (z<=y && argv[z][0] && argv[z][0]!='-')
3645 strncpy(sdf_path,argv[z],253);
3648 if (strcmp(argv[x],"-t")==0)
3650 /* Read Transmitter Location */
3654 while (z<=y && argv[z][0] && argv[z][0]!='-' && txsites<4)
3656 strncpy(txfile,argv[z],253);
3657 tx_site[txsites]=LoadQTH(txfile);
3664 if (strcmp(argv[x],"-L")==0)
3668 if (z<=y && argv[z][0] && argv[z][0]!='-')
3670 sscanf(argv[z],"%lf",&altitudeLR);
3673 fprintf(stdout,"c and L are exclusive options, ignoring L.\n");
3682 if (strcmp(argv[x],"-l")==0)
3686 if (z<=y && argv[z][0] && argv[z][0]!='-')
3688 strncpy(longley_file,argv[z],253);
3690 /* Doing this twice is harmless */
3695 if (strcmp(argv[x],"-r")==0)
3697 /* Read Receiver Location */
3701 if (z<=y && argv[z][0] && argv[z][0]!='-')
3703 strncpy(rxfile,argv[z],253);
3704 rx_site=LoadQTH(rxfile);
3708 if (strcmp(argv[x],"-s")==0)
3710 /* Read city file(s) */
3714 while (z<=y && argv[z][0] && argv[z][0]!='-' && cities<5)
3716 strncpy(city_file[cities],argv[z],253);
3723 if (strcmp(argv[x],"-b")==0)
3725 /* Read Boundary File(s) */
3729 while (z<=y && argv[z][0] && argv[z][0]!='-' && bfs<5)
3731 strncpy(boundary_file[bfs],argv[z],253);
3739 /* Perform some error checking on the arguments
3740 and switches parsed from the command-line.
3741 If an error is encountered, print a message
3742 and exit gracefully. */
3746 fprintf(stderr,"\n%c*** ERROR: No transmitter site(s) specified!\n\n",7);
3750 for (x=0, y=0; x<txsites; x++)
3752 if (tx_site[x].lat==0.0 && tx_site[x].lon==0.0)
3754 fprintf(stderr,"\n*** ERROR: Transmitter site #%d not found!",x+1);
3761 fprintf(stderr,"%c\n\n",7);
3765 if ((coverage+LRmap)==0 && rx_site.lat==0.0 && rx_site.lon==0.0)
3767 fprintf(stderr,"\n%c*** ERROR: No receiver site found or specified!\n\n",7);
3771 /* No errors were detected. Whew! :-) */
3773 /* If no SDF path was specified on the command line (-d), check
3774 for a path specified in the $HOME/.splat_path file. If the
3775 file is not found, then sdf_path[] remains NULL, and the
3776 current working directory is assumed to contain the SDF
3781 sprintf(string,"%s/.splat_path",getenv("HOME"));
3782 fd=fopen(string,"r");
3786 fgets(string,253,fd);
3788 /* Remove <CR> and/or <LF> from string */
3790 for (x=0; string[x]!=13 && string[x]!=10 && string[x]!=0 && x<253; x++);
3793 strncpy(sdf_path,string,253);
3799 /* Ensure a trailing '/' is present in sdf_path */
3805 if (sdf_path[x-1]!='/' && x!=0)
3812 fprintf(stdout,"%s",header);
3823 rxlat=(unsigned char)floor(rx_site.lat);
3824 rxlon=(unsigned char)floor(rx_site.lon);
3831 else if (rxlat<min_lat)
3840 else if (rxlon<min_lon)
3850 for (y=0, z=0; z<txsites; z++)
3852 txlat=(unsigned char)floor(tx_site[z].lat);
3853 txlon=(unsigned char)floor(tx_site[z].lon);
3860 else if (txlat<min_lat)
3869 else if (txlon<min_lon)
3880 if (min_lat!=0 && min_lon!=0 && max_lat!=0 && max_lon!=0)
3882 for (y=min_lon; y<=max_lon; y++)
3883 for (x=min_lat; x<=max_lat; x++)
3885 sprintf(string,"%u:%u:%u:%u",x, x+1, y, y+1);
3892 for (z=0; z<txsites; z++)
3894 /* "Ball park" estimates used to load any additional
3895 SDF files required to conduct this analysis. */
3897 tx_range=sqrt(1.5*(tx_site[z].alt+GetElevation(tx_site[z])));
3898 rx_range=sqrt(1.5*altitude);
3900 /* deg_range determines the maximum
3901 amount of topo data we read */
3903 deg_range=(tx_range+rx_range)/69.0;
3905 /* max_range sets the maximum size of the
3906 analysis. A small, non-zero amount can
3907 be used to shrink the size of the analysis
3908 and limit the amount of topo data read by
3909 SPLAT! A very large number will only increase
3910 the width of the analysis, not the size of
3914 max_range=tx_range+rx_range;
3916 if (max_range<(tx_range+rx_range))
3917 deg_range=max_range/69.0;
3919 /* Prevent the demand for a really wide coverage
3920 from allocating more slots than are available
3925 case 2: deg_limit=0.25;
3928 case 4: deg_limit=0.5;
3931 case 9: deg_limit=1.0;
3934 case 16: deg_limit=2.0;
3937 case 25: deg_limit=3.0;
3940 if (tx_site[z].lat<70.0)
3941 deg_range_lon=deg_range/cos(deg2rad*tx_site[z].lat);
3943 deg_range_lon=deg_range/cos(deg2rad*70.0);
3945 /* Correct for squares in degrees not being square in miles */
3947 if (deg_range>deg_limit)
3948 deg_range=deg_limit;
3950 if (deg_range_lon>deg_limit)
3951 deg_range_lon=deg_limit;
3953 north_min=(unsigned char)floor(tx_site[z].lat-deg_range);
3954 north_max=(unsigned char)floor(tx_site[z].lat+deg_range);
3955 west_min=(unsigned char)floor(tx_site[z].lon-deg_range_lon);
3956 west_max=(unsigned char)floor(tx_site[z].lon+deg_range_lon);
3961 else if (north_min<min_lat)
3967 else if (west_min<min_lon)
3970 if (north_max>max_lat)
3973 if (west_max>max_lon)
3977 if (min_lat!=0 && min_lon!=0 && max_lat!=0 && max_lon!=0)
3979 for (y=min_lon; y<=max_lon; y++)
3980 for (x=min_lat; x<=max_lat; x++)
3982 sprintf(string,"%u:%u:%u:%u",x, x+1, y, y+1);
3990 /* "Ball park" estimates used to load any additional
3991 SDF files required to conduct this analysis. */
3993 tx_range=sqrt(1.5*(tx_site[0].alt+GetElevation(tx_site[0])));
3994 rx_range=sqrt(1.5*altitudeLR);
3997 tx_range=sqrt(5.0*tx_site[0].alt);
3998 rx_range=sqrt(5.0*altitudeLR);
4001 /* deg_range determines the maximum
4002 amount of topo data we read */
4004 deg_range=(tx_range+rx_range)/69.0;
4006 /* max_range sets the maximum size of the
4007 analysis. A small, non-zero amount can
4008 be used to shrink the size of the analysis
4009 and limit the amount of topo data read by
4010 SPLAT! A very large number will only increase
4011 the width of the analysis, not the size of
4015 max_range=tx_range+rx_range;
4017 if (max_range<(tx_range+rx_range))
4018 deg_range=max_range/69.0;
4020 /* Prevent the demand for a really wide coverage
4021 from allocating more slots than are available
4026 case 2: deg_limit=0.25;
4029 case 4: deg_limit=0.5;
4032 case 9: deg_limit=1.0;
4035 case 16: deg_limit=2.0;
4038 case 25: deg_limit=3.0;
4041 if (tx_site[0].lat<70.0)
4042 deg_range_lon=deg_range/cos(deg2rad*tx_site[0].lat);
4044 deg_range_lon=deg_range/cos(deg2rad*70.0);
4046 /* Correct for squares in degrees not being square in miles */
4048 if (deg_range>deg_limit)
4049 deg_range=deg_limit;
4051 if (deg_range_lon>deg_limit)
4052 deg_range_lon=deg_limit;
4054 north_min=(unsigned char)floor(tx_site[0].lat-deg_range);
4055 north_max=(unsigned char)floor(tx_site[0].lat+deg_range);
4056 west_min=(unsigned char)floor(tx_site[0].lon-deg_range_lon);
4057 west_max=(unsigned char)floor(tx_site[0].lon+deg_range_lon);
4062 else if (north_min<min_lat)
4068 else if (west_min<min_lon)
4071 if (north_max>max_lat)
4074 if (west_max>max_lon)
4077 if (min_lat!=0 && min_lon!=0 && max_lat!=0 && max_lon!=0)
4079 for (y=min_lon; y<=max_lon; y++)
4080 for (x=min_lat; x<=max_lat; x++)
4082 sprintf(string,"%u:%u:%u:%u",x, x+1, y, y+1);
4093 for (x=0; x<txsites; x++)
4095 PlotCoverage(tx_site[x],altitude);
4096 PlaceMarker(tx_site[x]);
4099 SiteReport(tx_site[x]);
4107 PlotLRMap(tx_site[0],altitudeLR);
4108 PlaceMarker(tx_site[0]);
4111 SiteReport(tx_site[0]);
4118 PlaceMarker(rx_site);
4120 for (x=0; x<txsites; x++)
4122 PlaceMarker(tx_site[x]);
4129 PlotPath(tx_site[x],rx_site,1);
4133 PlotPath(tx_site[x],rx_site,8);
4137 PlotPath(tx_site[x],rx_site,16);
4141 PlotPath(tx_site[x],rx_site,32);
4146 ObstructionReport(tx_site[x],rx_site,report);
4154 for (x=0; x<bfs; x++)
4155 LoadBoundaries(boundary_file[x]);
4160 for (x=0; x<cities; x++)
4161 LoadCities(city_file[x]);
4167 WritePPMLR(mapfile);
4174 for (x=0; terrain_file[x]!='.' && terrain_file[x]!=0 && x<80; x++);
4176 if (terrain_file[x]=='.') /* extension */
4179 for (y=1, z=x, x++; terrain_file[x]!=0 && x<253 && y<14; x++, y++)
4180 ext[y]=terrain_file[x];
4188 ext[0]=0; /* No extension */
4192 for (count=0; count<txsites; count++)
4194 sprintf(string,"%s-%c%s%c",terrain_file,'1'+count,ext,0);
4195 GraphTerrain(tx_site[count],rx_site,string);
4200 GraphTerrain(tx_site[0],rx_site,terrain_file);
4207 for (x=0; elevation_file[x]!='.' && elevation_file[x]!=0 && x<80; x++);
4209 if (elevation_file[x]=='.') /* extension */
4212 for (y=1, z=x, x++; elevation_file[x]!=0 && x<253 && y<14; x++, y++)
4213 ext[y]=elevation_file[x];
4216 elevation_file[z]=0;
4221 ext[0]=0; /* No extension */
4222 elevation_file[x]=0;
4225 for (count=0; count<txsites; count++)
4227 sprintf(string,"%s-%c%s%c",elevation_file,'1'+count,ext,0);
4228 GraphElevation(tx_site[count],rx_site,string);
4233 GraphElevation(tx_site[0],rx_site,elevation_file);
4240 for (x=0; height_file[x]!='.' && height_file[x]!=0 && x<80; x++);
4242 if (height_file[x]=='.') /* extension */
4245 for (y=1, z=x, x++; height_file[x]!=0 && x<253 && y<14; x++, y++)
4246 ext[y]=height_file[x];
4254 ext[0]=0; /* No extension */
4258 for (count=0; count<txsites; count++)
4260 sprintf(string,"%s-%c%s%c",height_file,'1'+count,ext,0);
4261 GraphHeight(tx_site[count],rx_site,string);
4266 GraphHeight(tx_site[0],rx_site,height_file);
4273 for (x=0; longley_file[x]!='.' && longley_file[x]!=0 && x<80; x++);
4275 if (longley_file[x]=='.') /* extension */
4278 for (y=1, z=x, x++; longley_file[x]!=0 && x<253 && y<14; x++, y++)
4279 ext[y]=longley_file[x];
4287 ext[0]=0; /* No extension */
4291 for (count=0; count<txsites; count++)
4293 sprintf(string,"%s-%c%s%c",longley_file,'1'+count,ext,0);
4294 GraphLongley(tx_site[count],rx_site,string);
4299 GraphLongley(tx_site[0],rx_site,longley_file);