X-Git-Url: https://git.gag.com/?p=debian%2Fsplat;a=blobdiff_plain;f=splat.cpp;fp=splat.cpp;h=4c801c02913af133d86f894db828b8584296056f;hp=3e17692aafa79b69d410cb869cc403462873f38d;hb=866d49fe6fd5051b29c0fcfc5d8e4f338fdfbe47;hpb=cae76b32deb53ddbfb94b44de132a72435f56e88 diff --git a/splat.cpp b/splat.cpp index 3e17692..4c801c0 100644 --- a/splat.cpp +++ b/splat.cpp @@ -1,6 +1,6 @@ /**************************************************************************** * SPLAT: An RF Signal Propagation Loss and Terrain Analysis Tool * -* Last update: 31-Mar-2006 * +* Last update: 21-Dec-2006 * ***************************************************************************** * Project started in 1997 by John A. Magliacane, KD2BD * ***************************************************************************** @@ -56,48 +56,52 @@ #define ARRAYSIZE 30025 #endif -char string[255], sdf_path[255], opened=0, *splat_version={"1.1.1"}; +char string[255], sdf_path[255], opened=0, *splat_version={"1.2.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; + METERS_PER_FOOT=0.3048, KM_PER_MILE=1.609344, earthradius, + max_range=0.0; int min_north=90, max_north=-90, min_west=360, max_west=-1, max_elevation=-32768, min_elevation=32768, bzerror, maxdB=230; -struct site { double lat; - double lon; - double alt; - char name[50]; - } site; - -struct path { float lat[ARRAYSIZE]; - float lon[ARRAYSIZE]; - float elevation[ARRAYSIZE]; - float distance[ARRAYSIZE]; - int length; - } path; - -struct dem { 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 LR { double eps_dielect; - double sgm_conductivity; - double eno_ns_surfref; - double frq_mhz; - double conf; - double rel; - int radio_climate; - int pol; - } LR; +unsigned char got_elevation_pattern=0, got_azimuth_pattern=0, metric=0; + +struct site { double lat; + double lon; + float alt; + char name[50]; + } site; + +struct path { double lat[ARRAYSIZE]; + double lon[ARRAYSIZE]; + double elevation[ARRAYSIZE]; + double distance[ARRAYSIZE]; + int length; + } path; + +struct dem { 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 LR { double eps_dielect; + double sgm_conductivity; + double eno_ns_surfref; + double frq_mhz; + double conf; + double rel; + int radio_climate; + int pol; + float antenna_pattern[361][1001]; + } LR; double elev_l[ARRAYSIZE+10]; @@ -138,12 +142,19 @@ 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; + + char sign; + int degrees, minutes, seconds; + double a, b, c, d; if (decimal<0.0) + { decimal=-decimal; + sign=-1; + } + + else + sign=1; a=floor(decimal); b=60.0*(decimal-a); @@ -161,7 +172,7 @@ char *dec2dms(double decimal) seconds=59; string[0]=0; - sprintf(string,"%d%c %d\' %d\"", degrees, 176, minutes, seconds); + sprintf(string,"%d%c %d\' %d\"", degrees*sign, 176, minutes, seconds); return (string); } @@ -173,8 +184,8 @@ int OrMask(double lat, double lon, int value) the mask based on the latitude and longitude of the area pointed to. */ - int x, y, indx, minlat, minlon; - char found; + int x, y, indx, minlat, minlon; + char found; minlat=(int)floor(lat); minlon=(int)floor(lon); @@ -213,9 +224,9 @@ double GetElevation(struct site 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; + char found; + int x, y, indx, minlat, minlon; + double elevation; elevation=-5000.0; @@ -238,12 +249,41 @@ double GetElevation(struct site location) return elevation; } +int AddElevation(double lat, double lon, double height) +{ + /* This function adds a user-defined terrain feature + (in meters AGL) to the digital elevation model data + in memory. Does nothing and returns 0 for locations + not found in memory. */ + + char found; + int x, y, indx, minlat, minlon; + + minlat=(int)floor(lat); + minlon=(int)floor(lon); + + x=(int)(1199.0*(lat-floor(lat))); + y=(int)(1199.0*(lon-floor(lon))); + + for (indx=0, found=0; indxcos_test_angle) + { + block=1; + first_obstruction_angle=((acos(cos_test_angle))/deg2rad)-90.0; + } + } + + if (block) + elevation=first_obstruction_angle; + + else + elevation=((acos(cos_xmtr_angle))/deg2rad)-90.0; + + path=temp; + + return elevation; +} + double AverageTerrain(struct site source, double azimuthx, double start_distance, double end_distance) { /* This function returns the average terrain calculated in @@ -418,9 +540,9 @@ double AverageTerrain(struct site source, double azimuthx, double start_distance 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; + 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; @@ -509,9 +631,9 @@ double haat(struct site antenna) 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; + 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, @@ -586,7 +708,6 @@ void PlaceMarker(struct site location) lat=location.lat; lon=location.lon; - if (latxmin && (LonDiff(lon,ymax)<0.0) && (LonDiff(lon,ymin)>0.0)) { p1=1.0/1200.0; @@ -788,9 +909,9 @@ double ReadBearing(char *input) embedded within the numbers expressed in the input string. Decimal seconds are permitted. */ - double seconds, bearing=0.0; - char string[20]; - int a, b, length, degrees, minutes; + double seconds, bearing=0.0; + char string[20]; + int a, b, length, degrees, minutes; /* Copy "input" to "string", and ignore any extra spaces that might be present in the process. */ @@ -800,7 +921,7 @@ double ReadBearing(char *input) for (a=0, b=0; a=0 && x<=360 && fd!=NULL) + { + azimuth[x]+=amplitude; + read_count[x]++; + } + + fgets(string,254,fd); + pointer=strchr(string,';'); + + if (pointer!=NULL) + *pointer=0; + + sscanf(string,"%f %f",&az, &litude); + + } while (feof(fd)==0); + + fclose(fd); + + + /* Handle 0=360 degree ambiguity */ + + if ((read_count[0]==0) && (read_count[360]!=0)) + { + read_count[0]=read_count[360]; + azimuth[0]=azimuth[360]; + } + + if ((read_count[0]!=0) && (read_count[360]==0)) + { + read_count[360]=read_count[0]; + azimuth[360]=azimuth[0]; + } + + /* Average pattern values in case more than + one was read for each degree of azimuth. */ + + for (x=0; x<=360; x++) + { + if (read_count[x]>1) + azimuth[x]/=(float)read_count[x]; + } + + /* Interpolate missing azimuths + to completely fill the array */ + + last_index=-1; + next_index=-1; + + for (x=0; x<=360; x++) + { + if (read_count[x]!=0) + { + if (last_index==-1) + last_index=x; + else + next_index=x; + } + + if (last_index!=-1 && next_index!=-1) + { + valid1=azimuth[last_index]; + valid2=azimuth[next_index]; + + span=next_index-last_index; + delta=(valid2-valid1)/(float)span; + + for (y=last_index+1; y=360) + y-=360; + + azimuth_pattern[y]=azimuth[x]; + } + + azimuth_pattern[360]=azimuth_pattern[0]; + + got_azimuth_pattern=255; + } + + /* Read and process .el file */ + + fd=fopen(elfile,"r"); + + if (fd!=NULL) + { + for (x=0; x<=10000; x++) + { + el_pattern[x]=0.0; + read_count[x]=0; + } + + /* Read mechanical tilt (degrees) and + tilt azimuth in degrees measured + clockwise from true North. */ + + fgets(string,254,fd); + pointer=strchr(string,';'); + + if (pointer!=NULL) + *pointer=0; + + sscanf(string,"%f %f",&mechanical_tilt, &tilt_azimuth); + + /* Read elevation (degrees) and corresponding + normalized field radiation pattern amplitude + (0.0 to 1.0) until EOF is reached. */ + + fgets(string,254,fd); + pointer=strchr(string,';'); + + if (pointer!=NULL) + *pointer=0; + + sscanf(string,"%f %f", &elevation, &litude); + + while (feof(fd)==0) + { + /* Read in normalized radiated field values + for every 0.01 degrees of elevation between + -10.0 and +90.0 degrees */ + + x=(int)rintf(100.0*(elevation+10.0)); + + if (x>=0 && x<=10000) + { + el_pattern[x]+=amplitude; + read_count[x]++; + } + + fgets(string,254,fd); + pointer=strchr(string,';'); + + if (pointer!=NULL) + *pointer=0; + + sscanf(string,"%f %f", &elevation, &litude); + } + + fclose(fd); + + /* Average the field values in case more than + one was read for each 0.01 degrees of elevation. */ + + for (x=0; x<=10000; x++) + { + if (read_count[x]>1) + el_pattern[x]/=(float)read_count[x]; + } + + /* Interpolate between missing elevations (if + any) to completely fill the array and provide + radiated field values for every 0.01 degrees of + elevation. */ + + last_index=-1; + next_index=-1; + + for (x=0; x<=10000; x++) + { + if (read_count[x]!=0) + { + if (last_index==-1) + last_index=x; + else + next_index=x; + } + + if (last_index!=-1 && next_index!=-1) + { + valid1=el_pattern[last_index]; + valid2=el_pattern[next_index]; + + span=next_index-last_index; + delta=(valid2-valid1)/(float)span; + + for (y=last_index+1; y=360) + y-=360; + + while (y<0) + y+=360; + + if (x<=180) + slant_angle[y]=-(tilt_increment*(90.0-xx)); + + if (x>180) + slant_angle[y]=-(tilt_increment*(xx-270.0)); + } + } + + slant_angle[360]=slant_angle[0]; /* 360 degree wrap-around */ + + for (w=0; w<=360; w++) + { + tilt=slant_angle[w]; + + /** Convert tilt angle to + an array index offset **/ + + y=(int)rintf(100.0*tilt); + + /* Copy shifted el_pattern[10001] field + values into elevation_pattern[361][1001] + at the corresponding azimuth, downsampling + (averaging) along the way in chunks of 10. */ + + for (x=y, z=0; z<=1000; x+=10, z++) + { + for (sum=0.0, a=0; a<10; a++) + { + b=a+x; + + if (b>=0 && b<=10000) + sum+=el_pattern[b]; + if (b<0) + sum+=el_pattern[0]; + if (b>10000) + sum+=el_pattern[10000]; + } + + elevation_pattern[w][z]=sum/10.0; + } + } + + got_elevation_pattern=255; + } + + for (x=0; x<=360; x++) + { + for (y=0; y<=1000; y++) + { + if (got_elevation_pattern) + elevation=elevation_pattern[x][y]; + else + elevation=1.0; + + if (got_azimuth_pattern) + az=azimuth_pattern[x]; + else + az=1.0; + + LR.antenna_pattern[x][y]=az*elevation; + } + } +} + int LoadSDF_SDF(char *name) { /* This function reads uncompressed SPLAT Data Files (.sdf) @@ -925,9 +1413,10 @@ int LoadSDF_SDF(char *name) 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; + 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]; @@ -1149,10 +1638,10 @@ int LoadSDF_BZ(char *name) 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; + 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]; @@ -1299,9 +1788,11 @@ int LoadSDF_BZ(char *name) fflush(stdout); return 1; } + else return -1; } + else return 0; } @@ -1317,9 +1808,9 @@ char LoadSDF(char *name) 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; + int x, y, indx, minlat, minlon, maxlat, maxlon; + char found, free_slot=0; + int return_value=-1; /* Try to load an uncompressed SDF first. */ @@ -1448,10 +1939,10 @@ void LoadCities(char *filename) 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; + int x, y, z; + char input[80], str[3][80]; + struct site city_site; + FILE *fd=NULL; fd=fopen(filename,"r"); @@ -1496,45 +1987,185 @@ void LoadCities(char *filename) fprintf(stdout,"Done!\n"); fflush(stdout); } + else fprintf(stderr,"*** ERROR: \"%s\": not found!\n",filename); } -void LoadBoundaries(char *filename) +void LoadUDT(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; + /* This function reads a file containing User-Defined Terrain + features for their addition to the digital elevation model + data used by SPLAT!. Elevations in the UDT file are evaluated + and then copied into a temporary file under /tmp. Then the + contents of the temp file are scanned, and if found to be unique, + are added to the ground elevations described by the digital + elevation data already loaded into memory. */ + + int i, x, y, z, fd=0; + char input[80], str[3][80], tempname[15], *pointer=NULL; + double latitude, longitude, height, templat, templon, + tempheight, one_pixel; + FILE *fd1=NULL, *fd2=NULL; + + strcpy(tempname,"/tmp/XXXXXX\0"); + one_pixel=1.0/1200.0; - fd=fopen(filename,"r"); + fd1=fopen(filename,"r"); - if (fd!=NULL) + if (fd1!=NULL) { - fgets(string,78,fd); + fd=mkstemp(tempname); + fd2=fopen(tempname,"w"); + + fgets(input,78,fd1); + + pointer=strchr(input,';'); + + if (pointer!=NULL) + *pointer=0; fprintf(stdout,"Reading \"%s\"... ",filename); fflush(stdout); - do + while (feof(fd1)==0) { - fgets(string,78,fd); - sscanf(string,"%lf %lf", &lon0, &lat0); - fgets(string,78,fd); + /* Parse line for latitude, longitude, height */ - do + for (x=0, y=0, z=0; x<78 && input[x]!=0 && z<3; x++) { - sscanf(string,"%lf %lf", &lon1, &lat1); - - lon0=fabs(lon0); - lon1=fabs(lon1); + if (input[x]!=',' && y<78) + { + str[z][y]=input[x]; + y++; + } + + else + { + str[z][y]=0; + z++; + y=0; + } + } + + latitude=ReadBearing(str[0]); + longitude=ReadBearing(str[1]); + + /* Remove and/or from antenna height string */ + + for (i=0; str[2][i]!=13 && str[2][i]!=10 && str[2][i]!=0; i++); + + str[2][i]=0; + + /* The terrain feature may be expressed in either + 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. Otherwise the height is interpreted + as being expressed in feet. */ + + for (i=0; str[2][i]!='M' && str[2][i]!='m' && str[2][i]!=0 && i<48; i++); + + if (str[2][i]=='M' || str[2][i]=='m') + { + str[2][i]=0; + height=rint(atof(str[2])); + } + + else + { + str[2][i]=0; + height=rint(3.28084*atof(str[2])); + } + + if (height>0.0) + fprintf(fd2,"%f, %f, %f\n",latitude, longitude, height); + + fgets(input,78,fd1); + + pointer=strchr(input,';'); + + if (pointer!=NULL) + *pointer=0; + } + + fclose(fd1); + fclose(fd2); + close(fd); + + fprintf(stdout,"Done!\n"); + fflush(stdout); + + fd1=fopen(tempname,"r"); + fd2=fopen(tempname,"r"); + + fscanf(fd1,"%lf, %lf, %lf", &latitude, &longitude, &height); + + for (y=0; feof(fd1)==0; y++) + { + rewind(fd2); + + fscanf(fd2,"%lf, %lf, %lf", &templat, &templon, &tempheight); + + for (x=0, z=0; feof(fd2)==0; x++) + { + if (x>y) + if (fabs(latitude-templat)<=one_pixel && fabs(longitude-templon)<=one_pixel) + z=1; + + fscanf(fd2,"%lf, %lf, %lf", &templat, &templon, &tempheight); + } + + if (z==0) + AddElevation(latitude, longitude, height); + + fscanf(fd1,"%lf, %lf, %lf", &latitude, &longitude, &height); + } + + fclose(fd1); + fclose(fd2); + unlink(tempname); + } + + 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; @@ -1575,12 +2206,12 @@ void ReadLRParm(char *txsite_filename) 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". */ + function to be 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; + double din; + char filename[255], string[80], *pointer=NULL; + int iin, ok=0, x; + FILE *fd=NULL, *outfile=NULL; /* Default parameters in case things go bad */ @@ -1605,20 +2236,6 @@ void ReadLRParm(char *txsite_filename) 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) @@ -1633,10 +2250,10 @@ void ReadLRParm(char *txsite_filename) { fgets(string,80,fd); - for (x=0; lookup[(int)string[x]] && x<20; x++) - string[x]=lookup[(int)string[x]]; + pointer=strchr(string,';'); - string[x]=0; + if (pointer!=NULL) + *pointer=0; ok=sscanf(string,"%lf", &din); @@ -1646,10 +2263,10 @@ void ReadLRParm(char *txsite_filename) fgets(string,80,fd); - for (x=0; lookup[(int)string[x]] && x<20; x++) - string[x]=lookup[(int)string[x]]; + pointer=strchr(string,';'); - string[x]=0; + if (pointer!=NULL) + *pointer=0; ok=sscanf(string,"%lf", &din); } @@ -1660,10 +2277,10 @@ void ReadLRParm(char *txsite_filename) fgets(string,80,fd); - for (x=0; lookup[(int)string[x]] && x<20; x++) - string[x]=lookup[(int)string[x]]; + pointer=strchr(string,';'); - string[x]=0; + if (pointer!=NULL) + *pointer=0; ok=sscanf(string,"%lf", &din); } @@ -1674,10 +2291,10 @@ void ReadLRParm(char *txsite_filename) fgets(string,80,fd); - for (x=0; lookup[(int)string[x]] && x<20; x++) - string[x]=lookup[(int)string[x]]; + pointer=strchr(string,';'); - string[x]=0; + if (pointer!=NULL) + *pointer=0; ok=sscanf(string,"%lf", &din); } @@ -1688,10 +2305,10 @@ void ReadLRParm(char *txsite_filename) fgets(string,80,fd); - for (x=0; lookup[(int)string[x]] && x<20; x++) - string[x]=lookup[(int)string[x]]; + pointer=strchr(string,';'); - string[x]=0; + if (pointer!=NULL) + *pointer=0; ok=sscanf(string,"%d", &iin); } @@ -1702,10 +2319,10 @@ void ReadLRParm(char *txsite_filename) fgets(string,80,fd); - for (x=0; lookup[(int)string[x]] && x<20; x++) - string[x]=lookup[(int)string[x]]; + pointer=strchr(string,';'); - string[x]=0; + if (pointer!=NULL) + *pointer=0; ok=sscanf(string,"%d", &iin); } @@ -1716,10 +2333,10 @@ void ReadLRParm(char *txsite_filename) fgets(string,80,fd); - for (x=0; lookup[(int)string[x]] && x<20; x++) - string[x]=lookup[(int)string[x]]; + pointer=strchr(string,';'); - string[x]=0; + if (pointer!=NULL) + *pointer=0; ok=sscanf(string,"%lf", &din); } @@ -1730,18 +2347,21 @@ void ReadLRParm(char *txsite_filename) fgets(string,80,fd); - for (x=0; lookup[(int)string[x]] && x<20; x++) - string[x]=lookup[(int)string[x]]; + pointer=strchr(string,';'); - string[x]=0; + if (pointer!=NULL) + *pointer=0; ok=sscanf(string,"%lf", &din); } + fclose(fd); + if (ok) + { LR.rel=din; - - fclose(fd); + LoadPAT(filename); + } } if (fd==NULL) @@ -1752,21 +2372,13 @@ void ReadLRParm(char *txsite_filename) 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); @@ -1778,98 +2390,6 @@ void ReadLRParm(char *txsite_filename) 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 @@ -1928,39 +2448,146 @@ void PlotPath(struct site source, struct site destination, char mask_value) } } -void PlotLRPath(struct site source, struct site destination) +void PlotLRPath(struct site source, struct site destination, FILE *fd) { - /* 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; + /* This function plots the RF path loss between source and + destination points based on the Longley-Rice propagation + model, taking into account antenna pattern data, if available. */ + + char block=0, strmode[100]; + int x, y, errnum; + double loss, azimuth, pattern=0.0, + source_alt, dest_alt, source_alt2, dest_alt2, + cos_xmtr_angle, cos_test_angle=0.0, test_alt, + elevation, distance=0.0, + four_thirds_earth; + struct site temp; ReadPath(source,destination); - elev_l[1]=0.04*METERS_PER_MILE; + + four_thirds_earth=EARTHRADIUS*(4.0/3.0); + + /* Copy elevations along path into the elev_l[] array. */ for (x=0; xcos_test_angle) + block=1; + } + + /* At this point, we have the elevation angle + to the first obstruction (if it exists). */ + } + + /* Determine attenuation for each point along the + path using Longley-Rice's point_to_point mode + starting at y=2 (number_of_points = 1), the + shortest distance terrain can play a role in + path loss. */ + + elev_l[0]=y-1; /* (number of points - 1) */ + + /* Distance between elevation samples */ + elev_l[1]=METERS_PER_MILE*(path.distance[y]-path.distance[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); + 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); + + if (block) + elevation=((acos(cos_test_angle))/deg2rad)-90.0; + + else + elevation=((acos(cos_xmtr_angle))/deg2rad)-90.0; + + temp.lat=path.lat[y]; + temp.lon=path.lon[y]; + + azimuth=(Azimuth(source,temp)); + + if (fd!=NULL) + { + /* Write path loss data to output file */ + + fprintf(fd,"%.7f, %.7f, %.3f, %.3f, %.2f\n",path.lat[y], path.lon[y], azimuth, elevation, loss); + } + + /* Integrate the antenna's radiation + pattern into the overall path loss. */ - /* Note: PASS BY REFERENCE ... loss and errnum are - passed by reference, and are only used in this - file by this function */ + x=(int)rint(10.0*(10.0-elevation)); + + if (x>=0 && x<=1000) + { + azimuth=rint(azimuth); + + pattern=(double)LR.antenna_pattern[(int)azimuth][x]; + + if (pattern!=0.0) + { + pattern=20.0*log10(pattern); + loss-=pattern; + } + } if (loss>225.0) loss=225.0; @@ -1975,7 +2602,7 @@ void PlotLRPath(struct site source, struct site destination) 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) + else if (GetMask(path.lat[y],path.lon[y])==0 && path.distance[y]>max_range) OrMask(path.lat[y],path.lon[y],1); } } @@ -2011,7 +2638,7 @@ void PlotCoverage(struct site source, double altitude) 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); + fprintf(stdout,"\nComputing line-of-sight coverage of %s with an RX antenna\nat %.2f %s AGL:\n\n 0%c to 25%c ",source.name,metric?altitude*METERS_PER_FOOT:altitude,metric?"meters":"feet",37,37); fflush(stdout); /* 18.75=1200 pixels/degree divided by 64 loops @@ -2151,7 +2778,7 @@ void PlotCoverage(struct site source, double altitude) } } -void PlotLRMap(struct site source, double altitude) +void PlotLRMap(struct site source, double altitude, char *plo_filename) { /* This function performs a 360 degree sweep around the transmitter site (source location), and plots the @@ -2166,6 +2793,7 @@ void PlotLRMap(struct site source, double altitude) struct site edge; float lat, lon, one_pixel; unsigned char symbol[4], x; + FILE *fd=NULL; one_pixel=1.0/1200.0; @@ -2177,9 +2805,20 @@ void PlotLRMap(struct site source, double altitude) 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); + + fprintf(stdout,"out to a radius\nof %.2f %s with an RX antenna at %.2f %s AGL:\n\n 0%c to 25%c ",metric?max_range*KM_PER_MILE:max_range,metric?"kilometers":"miles",metric?altitude*METERS_PER_FOOT:altitude,metric?"meters":"feet",37,37); fflush(stdout); + if (plo_filename[0]!=0) + fd=fopen(plo_filename,"wb"); + + if (fd!=NULL) + { + /* Write header information to output file */ + + fprintf(fd,"%d, %d\t; max_west, min_west\n%d, %d\t; max_north, min_north\n",max_west, min_west, max_north, min_north); + } + /* 18.75=1200 pixels/degree divided by 64 loops per progress indicator symbol (.oOo) printed. */ @@ -2194,7 +2833,7 @@ void PlotLRMap(struct site source, double altitude) edge.lon=lon; edge.alt=altitude; - PlotLRPath(source,edge); + PlotLRPath(source,edge,fd); count++; if (count==z) @@ -2222,7 +2861,7 @@ void PlotLRMap(struct site source, double altitude) edge.lon=min_west; edge.alt=altitude; - PlotLRPath(source,edge); + PlotLRPath(source,edge,fd); count++; if (count==z) @@ -2253,7 +2892,7 @@ void PlotLRMap(struct site source, double altitude) edge.lon=lon; edge.alt=altitude; - PlotLRPath(source,edge); + PlotLRPath(source,edge,fd); count++; if (count==z) @@ -2281,7 +2920,7 @@ void PlotLRMap(struct site source, double altitude) edge.lon=max_west; edge.alt=altitude; - PlotLRPath(source,edge); + PlotLRPath(source,edge,fd); count++; if (count==z) @@ -2297,11 +2936,14 @@ void PlotLRMap(struct site source, double altitude) } } + if (fd!=NULL) + fclose(fd); + fprintf(stdout,"\nDone!\n"); fflush(stdout); } -void WritePPM(char *filename) +void WritePPM(char *filename, unsigned char geo) { /* This function generates a topographic map in Portable Pix Map (PPM) format based on logarithmically scaled topology data, @@ -2310,7 +2952,7 @@ void WritePPM(char *filename) from its representation in dem[][] so that north points up and east points right in the image generated. */ - char mapfile[255]; + char mapfile[255], geofile[255]; unsigned char found, mask; unsigned width, height, output; int indx, x, y, x0=0, y0=0, minlat, minlon; @@ -2329,13 +2971,35 @@ void WritePPM(char *filename) else { for (x=0; filename[x]!='.' && filename[x]!=0 && x<250; x++) + { mapfile[x]=filename[x]; + geofile[x]=filename[x]; + } mapfile[x]='.'; + geofile[x]='.'; mapfile[x+1]='p'; + geofile[x+1]='g'; mapfile[x+2]='p'; + geofile[x+2]='e'; mapfile[x+3]='m'; + geofile[x+3]='o'; mapfile[x+4]=0; + geofile[x+4]=0; + } + + if (geo) + { + fd=fopen(geofile,"wb"); + + fprintf(fd,"FILENAME\t%s\n",mapfile); + fprintf(fd,"#\t\tX\tY\tLong\t\tLat\n"); + fprintf(fd,"TIEPOINT\t0\t0\t%d.000\t\t%d.000\n",(max_west<180?-max_west:360-max_west),max_north); + fprintf(fd,"TIEPOINT\t%u\t%u\t%d.000\t\t%d.000\n",width-1,height-1,(min_west<180?-min_west:360-min_west),min_north); + fprintf(fd,"IMAGESIZE\t%u\t%u\n",width,height); + fprintf(fd,"#\n# Auto Generated by SPLAT! v%s\n#\n",splat_version); + + fclose(fd); } fd=fopen(mapfile,"wb"); @@ -2482,7 +3146,7 @@ void WritePPM(char *filename) fflush(stdout); } -void WritePPMLR(char *filename) +void WritePPMLR(char *filename, unsigned char geo) { /* This function generates a topographic map in Portable Pix Map (PPM) format based on the content of flags held in the mask[][] @@ -2490,7 +3154,7 @@ void WritePPMLR(char *filename) 90 degrees from its representation in dem[][] so that north points up and east points right in the image generated. */ - char mapfile[255]; + char mapfile[255], geofile[255]; unsigned width, height, output; unsigned char found, mask, cityorcounty; int indx, x, y, t, t2, x0, y0, minlat, minlon, loss; @@ -2509,13 +3173,35 @@ void WritePPMLR(char *filename) else { for (x=0; filename[x]!='.' && filename[x]!=0 && x<250; x++) + { mapfile[x]=filename[x]; + geofile[x]=filename[x]; + } mapfile[x]='.'; + geofile[x]='.'; mapfile[x+1]='p'; + geofile[x+1]='g'; mapfile[x+2]='p'; + geofile[x+2]='e'; mapfile[x+3]='m'; + geofile[x+3]='o'; mapfile[x+4]=0; + geofile[x+4]=0; + } + + if (geo) + { + fd=fopen(geofile,"wb"); + + fprintf(fd,"FILENAME\t%s\n",mapfile); + fprintf(fd,"#\t\tX\tY\tLong\t\tLat\n"); + fprintf(fd,"TIEPOINT\t0\t0\t%d.000\t\t%d.000\n",(max_west<180?-max_west:360-max_west),max_north); + fprintf(fd,"TIEPOINT\t%u\t%u\t%d.000\t\t%.3f\n",width-1,height+29,(min_west<180?-min_west:360-min_west),(double)(min_north-0.025)); + fprintf(fd,"IMAGESIZE\t%u\t%u\n",width,height+30); + fprintf(fd,"#\n# Auto Generated by SPLAT! v%s\n#\n",splat_version); + + fclose(fd); } fd=fopen(mapfile,"wb"); @@ -2554,7 +3240,6 @@ void WritePPMLR(char *filename) { /* Text Labels - Black or Red */ - /* if ((mask&120) && (loss<=maxdB)) */ if ((mask&120) && (loss<=90)) fprintf(fd,"%c%c%c",0,0,0); else @@ -2797,18 +3482,24 @@ void GraphTerrain(struct site source, struct site destination, char *name) 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. */ + If no extension is found, .png is assumed. */ - int x, y, z; - char filename[255], term[15], ext[15]; - FILE *fd=NULL; + int x, y, z; + char filename[255], term[30], ext[15]; + FILE *fd=NULL; ReadPath(destination,source); fd=fopen("profile.gp","wb"); for (x=0; xmaxangle) 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); + if (metric) + { + fprintf(fd,"%f\t%f\n",KM_PER_MILE*path.distance[path.length-1],refangle); + fprintf(fd2,"%f\t%f\n",KM_PER_MILE*path.distance[path.length-1],refangle); + } + + else + { + 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); @@ -2932,8 +3656,8 @@ void GraphElevation(struct site source, struct site destination, char *name) /* Default filename and output file type */ strncpy(filename,"profile\0",8); - strncpy(term,"gif\0",4); - strncpy(ext,"gif\0",4); + strncpy(term,"png\0",4); + strncpy(ext,"png\0",4); } else @@ -2957,11 +3681,11 @@ void GraphElevation(struct site source, struct site destination, char *name) } else - { /* No extension -- Default is gif */ + { /* No extension -- Default is png */ filename[x]=0; - strncpy(term,"gif\0",4); - strncpy(ext,"gif\0",4); + strncpy(term,"png\0",4); + strncpy(ext,"png\0",4); } } @@ -2972,7 +3696,7 @@ void GraphElevation(struct site source, struct site destination, char *name) strncpy(ext,"ps\0",3); else if (strncmp(ext,"ps",2)==0) - strncpy(term,"postscript\0",11); + strncpy(term,"postscript enhanced color\0",26); fprintf(stdout,"Writing \"%s.%s\"...",filename,ext); fflush(stdout); @@ -2981,12 +3705,19 @@ void GraphElevation(struct site source, struct site destination, char *name) 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 encoding iso_8859_1\n"); 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 title \"SPLAT! Elevation Profile Between %s and %s (%.2f%c azimuth)\"\n",destination.name,source.name,Azimuth(destination,source),176); + + if (metric) + fprintf(fd,"set xlabel \"Distance Between %s and %s (%.2f kilometers)\"\n",destination.name,source.name,KM_PER_MILE*Distance(source,destination)); + else + fprintf(fd,"set xlabel \"Distance Between %s and %s (%.2f miles)\"\n",destination.name,source.name,Distance(source,destination)); + + + fprintf(fd,"set ylabel \"Elevation Angle Along LOS 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"); + fprintf(fd,"plot \"profile.gp\" title \"Real Earth Profile\" with lines, \"reference.gp\" title \"Line of Sight Path (%.2f%c elevation)\" with lines\n",refangle,176); fclose(fd); @@ -3006,68 +3737,209 @@ void GraphElevation(struct site source, struct site destination, char *name) fprintf(stderr,"\n*** ERROR: Error occurred invoking gnuplot!\n"); } -void GraphHeight(struct site source, struct site destination, char *name) +void GraphHeight(struct site source, struct site destination, char *name, double f, unsigned char n) { /* 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; + and destination locations referenced to the line-of-sight path + between the receive 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, .png is assumed. */ + + int x, y, z; + char filename[255], term[30], ext[15]; + double a, b, c, height=0.0, refangle, cangle, maxheight=-100000.0, + minheight=100000.0, lambda=0.0, f_zone=0.0, fpt6_zone=0.0, + nm=0.0, nb=0.0, ed=0.0, es=0.0, r=0.0, d=0.0, d1=0.0, + terrain, azimuth, distance, dheight=0.0, minterrain=100000.0, + minearth=100000.0, miny, maxy, min2y, max2y; + struct site remote; + FILE *fd=NULL, *fd2=NULL, *fd3=NULL, *fd4=NULL, *fd5=NULL; ReadPath(destination,source); /* destination=RX, source=TX */ + azimuth=Azimuth(destination,source); + distance=Distance(destination,source); refangle=ElevationAngle(destination,source); b=GetElevation(destination)+destination.alt+earthradius; + /* Wavelength and path distance (great circle) in feet. */ + + if (f) + { + lambda=9.8425e8/(f*1e6); + d=5280.0*path.distance[path.length-1]; + } + + if (n) + { + ed=GetElevation(destination); + es=GetElevation(source); + nb=-destination.alt-ed; + nm=(-source.alt-es-nb)/(path.distance[path.length-1]); + } + fd=fopen("profile.gp","wb"); fd2=fopen("reference.gp","wb"); + fd5=fopen("curvature.gp", "wb"); - for (x=1; x0) + { + f_zone+=r; + fpt6_zone+=r; + } + } + + else + r=0.0; + + if (metric) + { + fprintf(fd,"%f\t%f\n",KM_PER_MILE*path.distance[x],METERS_PER_FOOT*height); + fprintf(fd2,"%f\t%f\n",KM_PER_MILE*path.distance[x],METERS_PER_FOOT*r); + fprintf(fd5,"%f\t%f\n",KM_PER_MILE*path.distance[x],METERS_PER_FOOT*(height-terrain)); + } + + else + { + fprintf(fd,"%f\t%f\n",path.distance[x],height); + fprintf(fd2,"%f\t%f\n",path.distance[x],r); + fprintf(fd5,"%f\t%f\n",path.distance[x],height-terrain); + } + + if (f) + { + if (metric) + { + fprintf(fd3,"%f\t%f\n",KM_PER_MILE*path.distance[x],METERS_PER_FOOT*f_zone); + fprintf(fd4,"%f\t%f\n",KM_PER_MILE*path.distance[x],METERS_PER_FOOT*fpt6_zone); + } + + else + { + fprintf(fd3,"%f\t%f\n",path.distance[x],f_zone); + fprintf(fd4,"%f\t%f\n",path.distance[x],fpt6_zone); + } + + if (f_zonemaxheight) maxheight=height; if (heightmaxheight) + maxheight=r; + + if (terrainmaxheight) + maxheight=r; - 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); + if (r-4999.0) - fprintf(fd2,"Antenna height above average terrain: %.2f feet\n",haavt); + { + if (metric) + fprintf(fd2,"Antenna height above average terrain: %.2f meters\n",METERS_PER_FOOT*haavt); + else + fprintf(fd2,"Antenna height above average terrain: %.2f feet\n",haavt); + } + + azimuth=Azimuth(source,destination); + angle1=ElevationAngle(source,destination); + angle2=ElevationAngle2(source,destination,earthradius); + + if (got_azimuth_pattern || got_elevation_pattern) + { + x=(int)rint(10.0*(10.0-angle2)); + + if (x>=0 && x<=1000) + pattern=(double)LR.antenna_pattern[(int)rint(azimuth)][x]; + + patterndB=20.0*log10(pattern); + + fprintf(fd2,"Antenna pattern between %s and %s: %.3f (%.2f dB)\n",source.name, destination.name, pattern, patterndB); + + } + + if (metric) + fprintf(fd2,"Distance to %s: %.2f kilometers\n",destination.name,METERS_PER_FOOT*Distance(source,destination)); - 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)); + else + fprintf(fd2,"Distance to %s: %.2f miles.\n",destination.name,Distance(source,destination)); + + fprintf(fd2,"Azimuth to %s: %.2f degrees\n",destination.name,azimuth); - angle=ElevationAngle(source,destination); + if (angle1>=0.0) + fprintf(fd2,"Elevation angle to %s: %+.4f degrees\n",destination.name,angle1); - if (angle>=0.0) - fprintf(fd2,"Angle of elevation between %s and %s: %+.4f degrees.\n",source.name,destination.name,angle); + else + fprintf(fd2,"Depression angle to %s: %+.4f degrees\n",destination.name,angle1); - if (angle<0.0) - fprintf(fd2,"Angle of depression between %s and %s: %+.4f degrees.\n",source.name,destination.name,angle); + if (angle1!=angle2) + { + if (angle2<0.0) + fprintf(fd2,"Depression"); + else + fprintf(fd2,"Elevation"); + + fprintf(fd2," angle to the first obstruction: %+.4f degrees\n",angle2); + } fprintf(fd2,"\n-------------------------------------------------------------------------\n\n"); @@ -3225,24 +4205,57 @@ void GraphLongley(struct site source, struct site destination, char *name) } 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)); + + if (metric) + { + fprintf(fd2,"Ground elevation: %.2f meters AMSL\n",METERS_PER_FOOT*GetElevation(destination)); + fprintf(fd2,"Antenna height: %.2f meters AGL / %.2f meters AMSL\n",METERS_PER_FOOT*destination.alt, METERS_PER_FOOT*(destination.alt+GetElevation(destination))); + } + + else + { + 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); + { + if (metric) + fprintf(fd2,"Antenna height above average terrain: %.2f meters\n",METERS_PER_FOOT*haavt); + else + fprintf(fd2,"Antenna height above average terrain: %.2f feet\n",haavt); + } + + if (metric) + fprintf(fd2,"Distance to %s: %.2f kilometers\n",source.name,KM_PER_MILE*Distance(source,destination)); + + else + fprintf(fd2,"Distance to %s: %.2f miles\n",source.name,Distance(source,destination)); + + azimuth=Azimuth(destination,source); - 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)); + angle1=ElevationAngle(destination,source); + angle2=ElevationAngle2(destination,source,earthradius); - angle=ElevationAngle(destination,source); + fprintf(fd2,"Azimuth to %s: %.2f degrees.\n",source.name,azimuth); - if (angle>=0.0) - fprintf(fd2,"Angle of elevation between %s and %s: %+.4f degrees.\n",destination.name,source.name,angle); + if (angle1>=0.0) + fprintf(fd2,"Elevation angle to %s: %+.4f degrees\n",source.name,angle1); - if (angle<0.0) - fprintf(fd2,"Angle of depression between %s and %s: %+.4f degrees.\n",destination.name,source.name,angle); + else + fprintf(fd2,"Depression angle to %s: %+.4f degrees\n",source.name,angle1); + + if (angle1!=angle2) + { + if (angle2<0.0) + fprintf(fd2,"Depression"); + else + fprintf(fd2,"Elevation"); + + fprintf(fd2," angle to the first obstruction: %+.4f degrees\n",angle2); + } fprintf(fd2,"\n-------------------------------------------------------------------------\n\n"); @@ -3297,54 +4310,145 @@ void GraphLongley(struct site source, struct site destination, char *name) 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 */ + ReadPath(source, destination); /* source=TX, destination=RX */ - elev_l[1]=0.04*METERS_PER_MILE; + /* Copy elevations along path into the elev_l[] array. */ - for (x=1; xcos_test_angle) + block=1; + } + + /* At this point, we have the elevation angle + to the first obstruction (if it exists). */ + } + + /* Determine path loss for each point along the + path using Longley-Rice's point_to_point mode + starting at x=2 (number_of_points = 1), the + shortest distance terrain can play a role in + path loss. */ + + elev_l[0]=y-1; /* (number of points - 1) */ + + /* Distance between elevation samples */ + elev_l[1]=METERS_PER_MILE*(path.distance[y]-path.distance[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); + 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); + + if (block) + elevation=((acos(cos_test_angle))/deg2rad)-90.0; + else + elevation=((acos(cos_xmtr_angle))/deg2rad)-90.0; + + /* Integrate the antenna's radiation + pattern into the overall path loss. */ + + x=(int)rint(10.0*(10.0-elevation)); + + if (x>=0 && x<=1000) + { + pattern=(double)LR.antenna_pattern[(int)azimuth][x]; + + if (pattern!=0.0) + patterndB=20.0*log10(pattern); + } - /* Note: PASS BY REFERENCE ... loss and errnum are - passed by reference, and are only used in this - file by this function */ + else + patterndB=0.0; + + total_loss=loss-patterndB; + if (metric) + { + fprintf(fd,"%f\t%f\n",KM_PER_MILE*(path.distance[path.length-1]-path.distance[y]),total_loss); + fprintf(fd2,"%7.2f\t\t%7.2f\t\t %d\t%s\n",KM_PER_MILE*path.distance[y],total_loss, errnum, strmode); + } - 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); + else + { + fprintf(fd,"%f\t%f\n",path.distance[path.length-1]-path.distance[y],total_loss); + fprintf(fd2,"%7.2f\t\t%7.2f\t\t %d\t%s\n",path.distance[y],total_loss, errnum, strmode); + } errflag|=errnum; - - if (loss>maxloss) - maxloss=loss; - if (lossmaxloss) + maxloss=total_loss; + + if (total_loss-4999.0) - fprintf(fd,"Antenna height above average terrain: %.2f feet\n",haavt); + { + if (metric) + fprintf(fd,"Antenna height above average terrain: %.2f meters\n",METERS_PER_FOOT*haavt); + else + fprintf(fd,"Antenna height above average terrain: %.2f feet\n",haavt); + } + + pattern=1.0; + patterndB=0.0; + distance=Distance(xmtr,rcvr); + azimuth=Azimuth(xmtr,rcvr); + angle1=ElevationAngle(xmtr,rcvr); + angle2=ElevationAngle2(xmtr,rcvr,earthradius); + + if (got_azimuth_pattern || got_elevation_pattern) + { + x=(int)rint(10.0*(10.0-angle2)); + + if (x>=0 && x<=1000) + pattern=(double)LR.antenna_pattern[(int)rint(azimuth)][x]; + + if (pattern!=1.0) + { + fprintf(fd,"Antenna pattern toward %s: %.3f",rcvr.name,pattern); + patterndB=20.0*log10(pattern); + fprintf(fd," (%.2f dB)\n",patterndB); + } + } - 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)); + if (metric) + fprintf(fd,"Distance to %s: %.2f kilometers\n",rcvr.name,KM_PER_MILE*distance); + + else + fprintf(fd,"Distance to %s: %.2f miles\n",rcvr.name,distance); - angle=ElevationAngle(xmtr,rcvr); + fprintf(fd,"Azimuth to %s: %.2f degrees\n",rcvr.name,azimuth); + + if (angle1>=0.0) + fprintf(fd,"Elevation angle to %s: %+.4f degrees\n",rcvr.name,angle1); + + else + fprintf(fd,"Depression angle to %s: %+.4f degrees\n",rcvr.name,angle1); - if (angle>=0.0) - fprintf(fd,"Angle of elevation between %s and %s: %+.4f degrees.\n",xmtr.name,rcvr.name,angle); + if (angle1!=angle2) + { + if (angle2<0.0) + fprintf(fd,"Depression"); + else + fprintf(fd,"Elevation"); - if (angle<0.0) - fprintf(fd,"Angle of depression between %s and %s: %+.4f degrees.\n",xmtr.name,rcvr.name,angle); + fprintf(fd," angle to the first obstruction: %+.4f degrees\n",angle2); + } fprintf(fd,"\n-------------------------------------------------------------------------\n\n"); @@ -3513,85 +4684,465 @@ void ObstructionReport(struct site xmtr, struct site rcvr, char report) fprintf(fd, " (%s S / ", 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)); + fprintf(fd, "%s W)\n", dec2dms(rcvr.lon)); + + if (metric) + { + fprintf(fd,"Ground elevation: %.2f meters AMSL\n",METERS_PER_FOOT*GetElevation(rcvr)); + fprintf(fd,"Antenna height: %.2f meters AGL / %.2f meters AMSL\n",METERS_PER_FOOT*rcvr.alt, METERS_PER_FOOT*(rcvr.alt+GetElevation(rcvr))); + } + + else + { + 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) + { + if (metric) + fprintf(fd,"Antenna height above average terrain: %.2f meters\n",METERS_PER_FOOT*haavt); + else + fprintf(fd,"Antenna height above average terrain: %.2f feet\n",haavt); + } + + azimuth=Azimuth(rcvr,xmtr); + angle1=ElevationAngle(rcvr,xmtr); + angle2=ElevationAngle2(rcvr,xmtr,earthradius); + + if (metric) + fprintf(fd,"Distance to %s: %.2f kilometers\n",xmtr.name,KM_PER_MILE*distance); + else + fprintf(fd,"Distance to %s: %.2f miles\n",xmtr.name,distance); + + fprintf(fd,"Azimuth to %s: %.2f degrees\n",xmtr.name,azimuth); + + if (angle1>=0.0) + fprintf(fd,"Elevation to %s: %+.4f degrees\n",xmtr.name,angle1); + + else + fprintf(fd,"Depression angle to %s: %+.4f degrees\n",xmtr.name,angle1); + + if (angle1!=angle2) + { + if (angle2<0.0) + fprintf(fd,"Depression"); + else + fprintf(fd,"Elevation"); + + fprintf(fd," angle to the first obstruction: %+.4f degrees\n",angle2); + + } + + fprintf(fd,"\n-------------------------------------------------------------------------\n\n"); + + if (report=='y') + { + /* Generate profile of the terrain. Create the path + from transmitter to receiver because that's the + way the original los() function did it, and going + the other way can yield slightly different results. */ + + ReadPath(xmtr,rcvr); + h_r=GetElevation(rcvr)+rcvr.alt+earthradius; + h_r_f1=h_r; + h_r_fpt6=h_r; + h_r_orig=h_r; + h_t=GetElevation(xmtr)+xmtr.alt+earthradius; + d_tx=5280.0*Distance(rcvr,xmtr); + cos_tx_angle=((h_r*h_r)+(d_tx*d_tx)-(h_t*h_t))/(2.0*h_r*d_tx); + cos_tx_angle_f1=cos_tx_angle; + cos_tx_angle_fpt6=cos_tx_angle; + + if (f) + lambda=9.8425e8/(f*1e6); + + /* At each point along the path calculate the cosine + of a sort of "inverse elevation angle" at the receiver. + From the antenna, 0 deg. looks at the ground, and 90 deg. + is parallel to the ground. + + Start at the receiver. If this is the lowest antenna, + then terrain obstructions will be nearest to it. (Plus, + that's the way the original los() did it.) + + Calculate cosines only. That's sufficient to compare + angles and it saves the extra computational burden of + acos(). However, note the inverted comparison: if + acos(A) > acos(B), then B > A. */ + + for (x=path.length-1; x>0; x--) + { + site_x.lat=path.lat[x]; + site_x.lon=path.lon[x]; + site_x.alt=0.0; + + h_x=GetElevation(site_x)+earthradius; + d_x=5280.0*Distance(rcvr,site_x); + + /* Deal with the LOS path first. */ + + cos_test_angle=((h_r*h_r)+(d_x*d_x)-(h_x*h_x))/(2.0*h_r*d_x); + + if (cos_tx_angle>cos_test_angle) + { + if (h_r==h_r_orig) + fprintf(fd,"SPLAT! detected obstructions at:\n\n"); + + if (site_x.lat>=0.0) + { + if (metric) + fprintf(fd,"\t%.4f N, %.4f W, %5.2f kilometers, %6.2f meters AMSL\n",site_x.lat, site_x.lon, KM_PER_MILE*(d_x/5280.0), METERS_PER_FOOT*(h_x-earthradius)); + else + fprintf(fd,"\t%.4f N, %.4f W, %5.2f miles, %6.2f feet AMSL\n",site_x.lat, site_x.lon, d_x/5280.0, h_x-earthradius); + } + + else + { + if (metric) + fprintf(fd,"\t%.4f S, %.4f W, %5.2f kilometers, %6.2f meters AMSL\n",-site_x.lat, site_x.lon, KM_PER_MILE*(d_x/5280.0), METERS_PER_FOOT*(h_x-earthradius)); + else + + fprintf(fd,"\t%.4f S, %.4f W, %5.2f miles, %6.2f feet AMSL\n",-site_x.lat, site_x.lon, d_x/5280.0, h_x-earthradius); + } + } + + while (cos_tx_angle>cos_test_angle) + { + h_r+=1; + cos_test_angle=((h_r*h_r)+(d_x*d_x)-(h_x*h_x))/(2.0*h_r*d_x); + cos_tx_angle=((h_r*h_r)+(d_tx*d_tx)-(h_t*h_t))/(2.0*h_r*d_tx); + } + + if (f) + { + /* Now clear the first Fresnel zone, but don't + clutter the obstruction report. */ + + cos_tx_angle_f1=((h_r_f1*h_r_f1)+(d_tx*d_tx)-(h_t*h_t))/(2.0*h_r_f1*d_tx); + h_los=sqrt(h_r_f1*h_r_f1+d_x*d_x-2*h_r_f1*d_x*cos_tx_angle_f1); + h_f=h_los-sqrt(lambda*d_x*(d_tx-d_x)/d_tx); + + while (h_fh_r_orig) + { + if (metric) + sprintf(string,"\nAntenna at %s must be raised to at least %.2f meters AGL\nto clear all obstructions detected by SPLAT!\n",rcvr.name, METERS_PER_FOOT*(h_r-GetElevation(rcvr)-earthradius)); + else + sprintf(string,"\nAntenna at %s must be raised to at least %.2f feet AGL\nto clear all obstructions detected by SPLAT!\n",rcvr.name, h_r-GetElevation(rcvr)-earthradius); + } + + else + sprintf(string,"\nNo obstructions to LOS path due to terrain were detected by SPLAT!\n"); + + if (f) + { + if (h_r_fpt6>h_r_orig) + { + if (metric) + sprintf(string_fpt6,"\nAntenna at %s must be raised to at least %.2f meters AGL\nto clear 60%c of the first Fresnel zone.\n",rcvr.name, METERS_PER_FOOT*(h_r_fpt6-GetElevation(rcvr)-earthradius),37); + + else + sprintf(string_fpt6,"\nAntenna at %s must be raised to at least %.2f feet AGL\nto clear 60%c of the first Fresnel zone.\n",rcvr.name, h_r_fpt6-GetElevation(rcvr)-earthradius,37); + } + + else + sprintf(string_fpt6,"\n60%c of the first Fresnel zone is clear.\n",37); + + if (h_r_f1>h_r_orig) + { + if (metric) + sprintf(string_f1,"\nAntenna at %s must be raised to at least %.2f meters AGL\nto clear the first Fresnel zone.\n",rcvr.name, METERS_PER_FOOT*(h_r_f1-GetElevation(rcvr)-earthradius)); + + else + sprintf(string_f1,"\nAntenna at %s must be raised to at least %.2f feet AGL\nto clear the first Fresnel zone.\n",rcvr.name, h_r_f1-GetElevation(rcvr)-earthradius); + + } + + else + sprintf(string_f1,"\nThe first Fresnel zone is clear.\n\n"); + } + } + + fprintf(fd,"%s",string); + + if (f) + { + fprintf(fd,"%s",string_f1); + fprintf(fd,"%s",string_fpt6); + } + + fclose(fd); + + /* Display report summary on terminal */ + + /* Line-of-sight status */ + + fprintf(stdout,"%s",string); + + if (f) + { + /* Fresnel zone status */ + + fprintf(stdout,"%s",string_f1); + fprintf(stdout,"%s",string_fpt6); + } + + fprintf(stdout, "\nObstruction report written to: \"%s\"\n",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"); + + if (xmtr.lat>=0.0) + { + fprintf(fd,"Site location: %.4f North / %.4f West",xmtr.lat, xmtr.lon); + fprintf(fd, " (%s N / ",dec2dms(xmtr.lat)); + } + + else + { + fprintf(fd,"Site location: %.4f South / %.4f West",-xmtr.lat, xmtr.lon); + fprintf(fd, " (%s S / ",dec2dms(xmtr.lat)); + } + + fprintf(fd, "%s W)\n",dec2dms(xmtr.lon)); + + if (metric) + { + fprintf(fd,"Ground elevation: %.2f meters AMSL\n",METERS_PER_FOOT*GetElevation(xmtr)); + fprintf(fd,"Antenna height: %.2f meters AGL / %.2f meters AMSL\n",METERS_PER_FOOT*xmtr.alt, METERS_PER_FOOT*(xmtr.alt+GetElevation(xmtr))); + } + + else + { + 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) + { + if (metric) + fprintf(fd,"Antenna height above average terrain: %.2f meters\n\n",METERS_PER_FOOT*terrain); + else + 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) + { + if (metric) + fprintf(fd,"%.2f meters AMSL\n",METERS_PER_FOOT*terrain); + else + 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); +} + +void LoadTopoData(int max_lon, int min_lon, int max_lat, int min_lat) +{ + /* This function loads the SDF files required + to cover the limits of the region specified. */ + + int x, y, width, ymin, ymax; + + width=ReduceAngle(max_lon-min_lon); + + if ((max_lon-min_lon)<=180.0) + { + for (y=0; y<=width; y++) + for (x=min_lat; x<=max_lat; x++) + { + ymin=(int)(min_lon+(double)y); + + while (ymin<0) + ymin+=360; + + while (ymin>=360) + ymin-=360; + + ymax=ymin+1; + + while (ymax<0) + ymax+=360; + + while (ymax>=360) + ymax-=360; + + sprintf(string,"%d:%d:%d:%d",x, x+1, ymin, ymax); + LoadSDF(string); + } + } + + else + { + for (y=0; y<=width; y++) + for (x=min_lat; x<=max_lat; x++) + { + ymin=max_lon+y; + + while (ymin<0) + ymin+=360; + + while (ymin>=360) + ymin-=360; + + ymax=ymin+1; + + while (ymax<0) + ymax+=360; + + while (ymax>=360) + ymax-=360; + + sprintf(string,"%d:%d:%d:%d",x, x+1, ymin, ymax); + LoadSDF(string); + } + } +} + +int LoadPLI(char *filename) +{ + int error=0, max_west, min_west, max_north, min_north; + char string[80], *pointer=NULL; + double latitude=0.0, longitude=0.0, azimuth=0.0, elevation=0.0, + loss=0.0; + FILE *fd; + + fd=fopen(filename,"r"); - haavt=haat(rcvr); + if (fd!=NULL) + { + fgets(string,78,fd); + pointer=strchr(string,';'); - if (haavt>-4999.0) - fprintf(fd,"Antenna height above average terrain: %.2f feet\n",haavt); + if (pointer!=NULL) + *pointer=0; + + sscanf(string,"%d, %d",&max_west, &min_west); - 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)); + fgets(string,78,fd); + pointer=strchr(string,';'); - angle=ElevationAngle(rcvr,xmtr); + if (pointer!=NULL) + *pointer=0; - if (angle>=0.0) - fprintf(fd,"Angle of elevation between %s and %s: %+.4f degrees.\n",rcvr.name,xmtr.name,angle); + sscanf(string,"%d, %d",&max_north, &min_north); - if (angle<0.0) - fprintf(fd,"Angle of depression between %s and %s: %+.4f degrees.\n",rcvr.name,xmtr.name,angle); + fgets(string,78,fd); + pointer=strchr(string,';'); - fprintf(fd,"\n-------------------------------------------------------------------------\n\n"); + if (pointer!=NULL) + *pointer=0; - if (report=='y') - { - /* Write an Obstruction Report */ + LoadTopoData(max_west-1, min_west, max_north-1, min_north); - new_site=rcvr; - result=los(xmtr,rcvr); - result2=result; - result2.alt-=1; - block=result.name[0]; + fprintf(stdout,"\nReading \"%s\"... ",filename); + fflush(stdout); - if (block=='*') - fprintf(fd,"SPLAT! detected obstructions at:\n\n"); + fscanf(fd,"%lf, %lf, %lf, %lf, %lf",&latitude, &longitude, &azimuth, &elevation, &loss); - while (block=='*') + while (feof(fd)==0) { - if (result.lat!=result2.lat || result.lon!=result2.lon || result.alt!=result2.alt) - { - if (result.lat>=0.0) - fprintf(fd,"\t%.4f N, %.4f W, %5.2f miles, %6.2f feet AMSL.\n",result.lat, result.lon, Distance(rcvr,result), result.alt); - else + if (loss>225.0) + loss=225.0; - fprintf(fd,"\t%.4f S, %.4f W, %5.2f miles, %6.2f feet AMSL.\n",-result.lat, result.lon, Distance(rcvr,result), result.alt); - } + if (loss<75.0) + loss=75.0; - result2=result; - new_site.alt+=1.0; + loss-=75.0; + loss/=10.0; + loss+=1.0; - /* Can you hear me now? :-) */ + if (loss<=(double)maxdB) + OrMask(latitude,longitude,((unsigned char)(loss))<<3); - result=los(xmtr,new_site); - block=result.name[0]; + fscanf(fd,"%lf, %lf, %lf, %lf, %lf",&latitude, &longitude, &azimuth, &elevation, &loss); } - 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); - fclose(fd); + fprintf(stdout," Done!\n"); + fflush(stdout); + } - /* Display LOS status to terminal */ + else + error=1; - fprintf(stdout,"%sObstruction report written to: \"%s\"\n",string,report_name); - fflush(stdout); + return error; } -void SiteReport(struct site xmtr) +void WriteKML(struct site source, struct site destination) { - char report_name[80]; - double terrain; - int x, azi; - FILE *fd; + int x, y; + char block, report_name[80]; + double distance, rx_alt, tx_alt, cos_xmtr_angle, + azimuth, cos_test_angle, test_alt; + FILE *fd=NULL; - sprintf(report_name,"%s-site_report.txt",xmtr.name); + ReadPath(source,destination); + + sprintf(report_name,"%s-to-%s.kml",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) @@ -3599,99 +5150,260 @@ void SiteReport(struct site xmtr) 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"); + fprintf(fd,"\n"); + fprintf(fd,"\n",splat_version); + fprintf(fd,"\n"); + fprintf(fd,"SPLAT! Path\n"); + fprintf(fd,"1\n"); + fprintf(fd,"Path Between %s and %s\n",source.name,destination.name); - fprintf(fd,"---------------------------------------------------------------------------\n\n"); + fprintf(fd,"\n"); + fprintf(fd," %s\n",source.name); + fprintf(fd," \n"); + fprintf(fd," Transmit Site\n"); - if (xmtr.lat>=0.0) - { - fprintf(fd,"Site location: %.4f North / %.4f West",xmtr.lat, xmtr.lon); - fprintf(fd, " (%s N / ",dec2dms(xmtr.lat)); - } + if (source.lat>=0.0) + fprintf(fd,"
%s North
\n",dec2dms(source.lat)); + else + fprintf(fd,"
%s South
\n",dec2dms(source.lat)); + + fprintf(fd,"
%s West
\n",dec2dms(source.lon)); + azimuth=Azimuth(source,destination); + distance=Distance(source,destination); + + if (metric) + fprintf(fd,"
%.2f km",distance*KM_PER_MILE); else - { - fprintf(fd,"Site location: %.4f South / %.4f West",-xmtr.lat, xmtr.lon); - fprintf(fd, " (%s S / ",dec2dms(xmtr.lat)); - } + fprintf(fd,"
%.2f miles",distance); + + fprintf(fd," to %s
\n
toward an azimuth of %.2f%c
\n",destination.name,azimuth,176); + + fprintf(fd," \n"); + fprintf(fd," 1\n"); + fprintf(fd," \n"); + fprintf(fd," \n"); + fprintf(fd," 1\n"); + fprintf(fd," relativeToGround\n"); + fprintf(fd," %f,%f,30\n",(source.lon<180.0?-source.lon:360.0-source.lon),source.lat); + fprintf(fd," \n"); + fprintf(fd,"\n"); + + fprintf(fd,"\n"); + fprintf(fd," %s\n",destination.name); + fprintf(fd," \n"); + fprintf(fd," Receive Site\n"); - 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)); + if (destination.lat>=0.0) + fprintf(fd,"
%s North
\n",dec2dms(destination.lat)); + else + fprintf(fd,"
%s South
\n",dec2dms(destination.lat)); - terrain=haat(xmtr); + fprintf(fd,"
%s West
\n",dec2dms(destination.lon)); - if (terrain>-4999.0) + + if (metric) + fprintf(fd,"
%.2f km",distance*KM_PER_MILE); + else + fprintf(fd,"
%.2f miles",distance); + + fprintf(fd," to %s
\n
toward an azimuth of %.2f%c
\n",source.name,Azimuth(destination,source),176); + + fprintf(fd," \n"); + fprintf(fd," 1\n"); + fprintf(fd," \n"); + fprintf(fd," \n"); + fprintf(fd," 1\n"); + fprintf(fd," relativeToGround\n"); + fprintf(fd," %f,%f,30\n",(destination.lon<180.0?-destination.lon:360.0-destination.lon),destination.lat); + fprintf(fd," \n"); + fprintf(fd,"\n"); + + fprintf(fd,"\n"); + fprintf(fd,"Point-to-Point Path\n"); + fprintf(fd," 1\n"); + fprintf(fd," 0\n"); + fprintf(fd," \n"); + fprintf(fd," \n"); + fprintf(fd," 1\n"); + fprintf(fd," 1\n"); + fprintf(fd," relativeToGround\n"); + fprintf(fd," \n"); + + for (x=0; x\n"); + fprintf(fd," \n"); + fprintf(fd,"\n"); + + fprintf(fd,"\n"); + fprintf(fd,"Line-of-Sight Path\n"); + fprintf(fd," 1\n"); + fprintf(fd," 0\n"); + fprintf(fd," \n"); + fprintf(fd," \n"); + fprintf(fd," 1\n"); + fprintf(fd," 1\n"); + fprintf(fd," relativeToGround\n"); + fprintf(fd," \n"); + + /* Walk across the "path", indentifying obstructions along the way */ + + for (y=0; y=0 && block==0; x--) { - fprintf(fd,"Average terrain at %3d degrees azimuth: ",azi); - terrain=AverageTerrain(xmtr,(double)azi,2.0,10.0); + distance=5280.0*(path.distance[y]-path.distance[x]); + test_alt=earthradius+path.elevation[x]; - if (terrain>-4999.0) - fprintf(fd,"%.2f feet AMSL\n",terrain); - else - fprintf(fd,"No terrain\n"); + 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 + an obstruction 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) + fprintf(fd," %f,%f,-30\n",(path.lon[y]<180.0?-path.lon[y]:360.0-path.lon[y]),path.lat[y]); + else + fprintf(fd," %f,%f,5\n",(path.lon[y]<180.0?-path.lon[y]:360.0-path.lon[y]),path.lat[y]); } - fprintf(fd,"\n---------------------------------------------------------------------------\n\n"); + fprintf(fd," \n"); + fprintf(fd," \n"); + fprintf(fd,"\n"); + + fprintf(fd," \n"); + fprintf(fd," %f\n",(source.lon<180.0?-source.lon:360.0-source.lon)); + fprintf(fd," %f\n",source.lat); + fprintf(fd," 300.0\n"); + fprintf(fd," 45.0\n"); + fprintf(fd," %f\n",azimuth); + fprintf(fd," \n"); + + fprintf(fd,"\n"); + fprintf(fd,"\n"); + fclose(fd); - fprintf(stdout,"\nSite analysis report written to: \"%s\"\n",report_name); + + fprintf(stdout, "KML file written to: \"%s\"\n",report_name); + + fflush(stdout); } int main(char argc, char *argv[]) { - int x, y, ymin, ymax, width, z=0, min_lat, min_lon, - max_lat, max_lon, rxlat, rxlon, txlat, txlon, - west_min, west_max, north_min, north_max; + int x, y, z=0, min_lat, min_lon, max_lat, max_lon, + rxlat, rxlon, txlat, txlon, west_min, west_max, + north_min, north_max; unsigned char 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, report='y'; + count, report='y', norm=0, topomap=0, geo=0, + kml=0; 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], *env=NULL, txfile[255], map=0, boundary_file[5][255], - rxsite=0; + udt_file[255], rxsite=0, plo_filename[255], + pli_filename[255], nf=0; 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; + deg_range_lon, er_mult, freq=0.0; struct site tx_site[4], rx_site; FILE *fd; - sprintf(header,"\n\t\t--==[ Welcome To SPLAT! v%s ]==--\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"); - fprintf(stdout,"\t-db maximum loss contour to display on path loss maps (80-230 dB)\n\n"); - + fprintf(stdout,"\n\t\t --==[ SPLAT! v%s Available Options... ]==--\n\n",splat_version); + fprintf(stdout," -t txsite(s).qth (max of 4)\n"); + fprintf(stdout," -r rxsite.qth\n"); + fprintf(stdout," -c plot coverage of TX(s) with an RX antenna at X feet/meters AGL\n"); + fprintf(stdout," -L plot path loss map of TX based on an RX at X feet/meters AGL\n"); + fprintf(stdout," -s filename(s) of city/site file(s) to import (5 max)\n"); + fprintf(stdout," -b filename(s) of cartographic boundary file(s) to import (max of 5)\n"); + fprintf(stdout," -p filename of terrain profile graph to plot\n"); + fprintf(stdout," -e filename of terrain elevation graph to plot\n"); + fprintf(stdout," -h filename of terrain height graph to plot\n"); + fprintf(stdout," -H filename of normalized terrain height graph to plot\n"); + fprintf(stdout," -l filename of Longley-Rice graph to plot\n"); + fprintf(stdout," -o filename of topographic map to generate (.ppm)\n"); + fprintf(stdout," -u filename of user-defined terrain file to import\n"); + fprintf(stdout," -d sdf file directory path (overrides path in ~/.splat_path file)\n"); + fprintf(stdout," -n no analysis, brief report\n"); + fprintf(stdout," -N no analysis, no report\n"); + fprintf(stdout," -m earth radius multiplier\n"); + fprintf(stdout," -f frequency for Fresnel zone calculation (MHz)\n"); + fprintf(stdout," -R modify default range for -c or -L (miles/kilometers)\n"); + fprintf(stdout," -db maximum loss contour to display on path loss maps (80-230 dB)\n"); + fprintf(stdout," -nf do not plot Fresnel zones in height plots\n"); + fprintf(stdout," -plo filename of path-loss output file\n"); + fprintf(stdout," -pli filename of path-loss input file\n"); + fprintf(stdout," -udt filename of user defined terrain input file\n"); + fprintf(stdout," -geo generate an Xastir .geo georeference file (with .ppm output)\n"); + fprintf(stdout," -kml generate a Google Earth .kml file (for point-to-point links)\n"); + fprintf(stdout," -metric employ metric rather than imperial units for all user I/O\n\n"); + + fprintf(stdout,"If that flew by too fast, consider piping the output through 'less':\n"); + fprintf(stdout,"\n\tsplat | less\n\n"); fprintf(stdout,"Type 'man splat', or see the documentation for more details.\n\n"); fflush(stdout); return 1; @@ -3699,6 +5411,7 @@ int main(char argc, char *argv[]) y=argc-1; + metric=0; rxfile[0]=0; txfile[0]=0; string[0]=0; @@ -3706,11 +5419,17 @@ int main(char argc, char *argv[]) elevation_file[0]=0; terrain_file[0]=0; sdf_path[0]=0; + udt_file[0]=0; path.length=0; + LR.frq_mhz=0.0; rx_site.lat=91.0; rx_site.lon=361.0; + plo_filename[0]=0; + pli_filename[0]=0; earthradius=EARTHRADIUS; + sprintf(header,"\n\t\t--==[ Welcome To SPLAT! v%s ]==--\n\n", splat_version); + for (x=0; x<4; x++) { tx_site[x].lat=91.0; @@ -3727,7 +5446,6 @@ int main(char argc, char *argv[]) dem[x].max_west=-1; } - /* Scan for command line arguments */ for (x=1; x<=y; x++) @@ -3775,6 +5493,14 @@ int main(char argc, char *argv[]) map=1; } + if (strcmp(argv[x],"-u")==0) + { + z=x+1; + + if (z<=y && argv[z][0] && argv[z][0]!='-') + strncpy(udt_file,argv[z],253); + } + if (strcmp(argv[x],"-c")==0) { z=x+1; @@ -3786,7 +5512,7 @@ int main(char argc, char *argv[]) } } - if (strcmp(argv[x],"-db")==0) + if (strcmp(argv[x],"-db")==0 || strcmp(argv[x],"-dB")==0) { z=x+1; @@ -3826,7 +5552,7 @@ int main(char argc, char *argv[]) } } - if (strcmp(argv[x],"-h")==0) + if (strcmp(argv[x],"-h")==0 || strcmp(argv[x],"-H")==0) { z=x+1; @@ -3835,26 +5561,37 @@ int main(char argc, char *argv[]) strncpy(height_file,argv[z],253); height_plot=1; } + + if (strcmp(argv[x],"-H")==0) + norm=1; + else + norm=0; } if (strcmp(argv[x],"-n")==0) { - if (z<=y && argv[z][0] && argv[z][0]!='-') - { - report='n'; - map=1; - } + report='n'; + map=1; } if (strcmp(argv[x],"-N")==0) { - if (z<=y && argv[z][0] && argv[z][0]!='-'); - { - report='N'; - map=1; - } + report='N'; + map=1; } + if (strcmp(argv[x],"-metric")==0) + metric=1; + + if (strcmp(argv[x],"-geo")==0) + geo=1; + + if (strcmp(argv[x],"-kml")==0) + kml=1; + + if (strcmp(argv[x],"-nf")==0) + nf=1; + if (strcmp(argv[x],"-d")==0) { z=x+1; @@ -3876,6 +5613,7 @@ int main(char argc, char *argv[]) txsites++; z++; } + z--; } @@ -3936,6 +5674,7 @@ int main(char argc, char *argv[]) cities++; z++; } + z--; } @@ -3951,8 +5690,41 @@ int main(char argc, char *argv[]) bfs++; z++; } + z--; } + + if (strcmp(argv[x],"-f")==0) + { + z=x+1; + + if (z<=y && argv[z][0] && argv[z][0]!='-') + { + sscanf(argv[z],"%lf",&freq); + + if (freq<20) + freq=20; + + if (freq>20e3) + freq=20e3; + } + } + + if (strcmp(argv[x],"-plo")==0) + { + z=x+1; + + if (z<=y && argv[z][0] && argv[z][0]!='-') + strncpy(plo_filename,argv[z],253); + } + + if (strcmp(argv[x],"-pli")==0) + { + z=x+1; + + if (z<=y && argv[z][0] && argv[z][0]!='-') + strncpy(pli_filename,argv[z],253); + } } /* Perform some error checking on the arguments @@ -3981,13 +5753,34 @@ int main(char argc, char *argv[]) exit (-1); } - if ((coverage+LRmap)==0 && rx_site.lat==91.0 && rx_site.lon==361.0) + if ((coverage+LRmap+pli_filename[0])==0 && rx_site.lat==91.0 && rx_site.lon==361.0) { - fprintf(stderr,"\n%c*** ERROR: No receiver site found or specified!\n\n",7); - exit (-1); + if (max_range!=0.0 && txsites!=0) + { + /* Plot topographic map of radius "max_range" */ + + map=0; + topomap=1; + report='N'; + } + + else + { + fprintf(stderr,"\n%c*** ERROR: No receiver site found or specified!\n\n",7); + exit (-1); + } } - /* No errors were detected. Whew! :-) */ + /* No major errors were detected. Whew! :-) */ + + /* Adjust input parameters if -metric option is used */ + + if (metric) + { + altitudeLR/=METERS_PER_FOOT; /* meters --> feet */ + max_range/=KM_PER_MILE; /* kilometers --> miles */ + altitude/=METERS_PER_FOOT; /* meters --> feet */ + } /* If no SDF path was specified on the command line (-d), check for a path specified in the $HOME/.splat_path file. If the @@ -4032,6 +5825,33 @@ int main(char argc, char *argv[]) fprintf(stdout,"%s",header); fflush(stdout); + if (pli_filename[0]) + { + y=LoadPLI(pli_filename); + + for (x=0; x=360) - ymin-=360; - - ymax=ymin+1; - - while (ymax<0) - ymax+=360; - - while (ymax>=360) - ymax-=360; - - sprintf(string,"%d:%d:%d:%d",x, x+1, ymin, ymax); - LoadSDF(string); - } - } - - else - { - for (y=0; y<=width; y++) - for (x=min_lat; x<=max_lat; x++) - { - ymin=max_lon+y; - - while (ymin<0) - ymin+=360; - - while (ymin>=360) - ymin-=360; - - ymax=ymin+1; - - while (ymax<0) - ymax+=360; - - while (ymax>=360) - ymax-=360; - - sprintf(string,"%d:%d:%d:%d",x, x+1, ymin, ymax); - LoadSDF(string); - } - } + LoadTopoData(max_lon, min_lon, max_lat, min_lat); - if (coverage | LRmap) + if (coverage | LRmap | topomap) { if (LRmap) txsites=1; @@ -4237,68 +6004,22 @@ int main(char argc, char *argv[]) max_lon=west_max; } + /* Load any additional SDF files, if required */ - /* Load the required SDF files */ - - width=ReduceAngle(max_lon-min_lon); - - if ((max_lon-min_lon)<=180.0) - { - for (y=0; y<=width; y++) - for (x=min_lat; x<=max_lat; x++) - { - ymin=(int)(min_lon+(double)y); - - while (ymin<0) - ymin+=360; - - while (ymin>=360) - ymin-=360; - - ymax=ymin+1; - - while (ymax<0) - ymax+=360; - - while (ymax>=360) - ymax-=360; - - sprintf(string,"%d:%d:%d:%d",x, x+1, ymin, ymax); - LoadSDF(string); - } - } - - else - { - for (y=0; y<=width; y++) - for (x=min_lat; x<=max_lat; x++) - { - ymin=(int)(max_lon+(double)y); - - while (ymin<0) - ymin+=360; - - while (ymin>=360) - ymin-=360; - - ymax=ymin+1; - - while (ymax<0) - ymax+=360; - - while (ymax>=360) - ymax-=360; - - sprintf(string,"%d:%d:%d:%d",x, x+1, ymin, ymax); - LoadSDF(string); - } - } + LoadTopoData(max_lon, min_lon, max_lat, min_lat); } + if (udt_file[0]) + LoadUDT(udt_file); - if (mapfile[0]) + if (mapfile[0] && topomap==0) map=1; + if (freq==0.0 && nf==0) + freq=LR.frq_mhz; + else + freq=0.0; + if (coverage | LRmap) { for (x=0; x1)