3.04 Egli model and unit test

CHANGELOG

@@ -1,5 +1,11 @@
 SIGNAL SERVER CHANGELOG 
 
+3.04 - 19 April 2017
+Added Egli VHF/UHF model courtesy of G6DTX
+Adjusted SUI correction factor for height. Most academic papers have /2000 but some are /2. Only out by 1e3 :O 
+New propagation model unit test script at models/test.cc
+Added 10 to ARRAYSIZE again as it was needed :p
+
 3.03 - 19 March 2017
 Path profile bugfix for some models
 Error handling for when prop loss < free space loss

Makefile

@@ -8,7 +8,7 @@ LIBS		= -lm -lpthread -ldl
 
 VPATH		= models
 objects 	= main.o cost.o ecc33.o ericsson.o fspl.o hata.o itwom3.0.o \
-		  los.o sui.o pel.o inputs.o outputs.o image.o image-ppm.o
+		  los.o sui.o pel.o egli.o inputs.o outputs.o image.o image-ppm.o
 
 GCC_MAJOR	:= $(shell $(CXX) -dumpversion 2>&1 | cut -d . -f 1)
 GCC_MINOR	:= $(shell $(CXX) -dumpversion 2>&1 | cut -d . -f 2)
@@ -43,14 +43,14 @@ main.o: main.cc common.h inputs.hh outputs.hh itwom3.0.hh los.hh
 inputs.o: inputs.cc common.h main.hh
 
 outputs.o: outputs.cc common.h inputs.hh main.hh cost.hh ecc33.hh ericsson.hh \
-	   fspl.hh hata.hh itwom3.0.hh sui.hh pel.hh
+	   fspl.hh hata.hh itwom3.0.hh sui.hh pel.hh egli.hh
 
 image.o: image.cc image-ppm.o
 
 image-ppm.o: image-ppm.cc
 
 los.o: los.cc common.h main.hh cost.hh ecc33.hh ericsson.hh fspl.hh hata.hh \
-       itwom3.0.hh sui.hh pel.hh
+       itwom3.0.hh sui.hh pel.hh egli.hh
 
 .PHONY: clean
 clean:

main.cc

@@ -1,9 +1,8 @@
-double version = 3.03;
+double version = 3.04;
 /****************************************************************************\
 *  Signal Server: Radio propagation simulator by Alex Farrant QCVS, 2E0TDW   *
 ******************************************************************************
 *    SPLAT! Project started in 1997 by John A. Magliacane, KD2BD             *
-*                                                                            *
 ******************************************************************************
 *         Please consult the SPLAT! documentation for a complete list of     *
 *         individuals who have contributed to this project.                  *
@@ -40,7 +39,7 @@ double version = 3.03;
 
 int MAXPAGES = 10*10;
 int IPPD = 1200;
-int ARRAYSIZE = MAXPAGES * IPPD;
+int ARRAYSIZE = (MAXPAGES * IPPD) + 10;
 
 char string[255], sdf_path[255], udt_file[255], opened = 0, gpsav =
     0, ss_name[16], dashes[80];
@@ -1105,7 +1104,7 @@ int main(int argc, char *argv[])
 		fprintf(stdout, "     -R Radius (miles/kilometers)\n");
 		fprintf(stdout,	"     -res Pixels per tile. 300/600/1200/3600 (Optional. LIDAR res is within the tile)\n");
 		fprintf(stdout,	"     -pm Propagation model. 1: ITM, 2: LOS, 3: Hata, 4: ECC33,\n");
-		fprintf(stdout,	"     	  5: SUI, 6: COST-Hata, 7: FSPL, 8: ITWOM, 9: Ericsson, 10: Plane earth\n");
+		fprintf(stdout,	"     	  5: SUI, 6: COST-Hata, 7: FSPL, 8: ITWOM, 9: Ericsson, 10: Plane earth, 11: Egli VHF/UHF\n");
 		fprintf(stdout,	"     -pe Propagation model mode: 1=Urban,2=Suburban,3=Rural\n");
 		fprintf(stdout,	"     -ked Knife edge diffraction (Already on for ITM)\n");
 		fprintf(stdout, "Debugging:\n");
@@ -1934,7 +1933,8 @@ int main(int argc, char *argv[])
 		PlotPath(tx_site[0], tx_site[1], 1);
 		PathReport(tx_site[0], tx_site[1], tx_site[0].filename, 0,
 			   propmodel, pmenv, rxGain);
-		SeriesData(tx_site[0], tx_site[1], tx_site[0].filename, 1,
+		// Order flipped for benefit of graph. Makes no difference to data.
+		SeriesData(tx_site[1], tx_site[0], tx_site[0].filename, 1,
 			   normalise);
 	}
 	fflush(stderr);

models/cost.cc

@@ -27,12 +27,12 @@ Distance 1-20km
 modes 1 = URBAN, 2 = SUBURBAN, 3 = OPEN
 http://morse.colorado.edu/~tlen5510/text/classwebch3.html
 */
-	if (f < 150 || f > 2000) {
+/*	if (f < 150 || f > 2000) {
 		fprintf
 		    (stderr,"Error: COST231 Hata model frequency range 150-2000MHz\n");
 		exit(EXIT_FAILURE);
 	}
-
+*/
 	int C = 3;		// 3dB for Urban
 	float lRxH = log10(11.75 * RxH);
 	float C_H = 3.2 * (lRxH * lRxH) - 4.97;	// Large city (conservative)

models/ecc33.cc

@@ -10,10 +10,11 @@ double ECC33pathLoss(float f, float TxH, float RxH, float d, int mode)
 		RxH=RxH/(d*2);
 	}
 
-	if (f < 700 || f > 3500) {
+/*	if (f < 700 || f > 3500) {
 		fprintf(stderr,"Error: ECC33 model frequency range 700-3500MHz\n");
 		exit(EXIT_FAILURE);
 	}
+*/
 	// MHz to GHz
 	f = f / 1000;
 

models/egli.cc

@@ -0,0 +1,83 @@
+/*****************************************************************************
+*  Egli VHF/UHF model for Signal Server by G6DTX                             *
+*  April 2017                                                                *
+*  This program is free software; you can redistribute it and/or modify it   *
+*  under the terms of the GNU General Public License as published by the     *
+*  Free Software Foundation; either version 2 of the License or any later    *
+*  version.                                                                  *
+*                                                                            *
+*  This program is distributed in the hope that it will useful, but WITHOUT  *
+*  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or     *
+*  FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License v3  *
+*  for more details.                                                         *
+******************************************************************************
+
+Frequency 30 to 1000MHz
+h1 = 1m and above
+h2 = 1m and above
+Distance 1 to 50km
+
+http://people.seas.harvard.edu/~jones/es151/prop_models/propagation.html#pel
+*/
+
+#include <stdio.h>
+#include <stdlib.h>
+#include <math.h>
+
+//static float fcmin = 30.0;
+//static float fcmax = 1000.0;
+//static float dmin  = 1.0;
+//static float dmax  = 50.0;
+//static float h1min = 1.0;
+//static float h2min = 1.0;
+
+static __inline float _10log10f(float x)
+{
+  return(4.342944f*logf(x));
+}
+
+double EgliPathLoss(float f, float h1, float h2, float d)
+{
+  double Lp50 = NAN;
+  float C1, C2;
+
+/*  if ((f >= fcmin) && (f <= fcmax) &&
+      (h1 >= h1min) && (h2 >= h2min))
+  {*/
+    if (h1 > 10.0 && h2 > 10.0)
+    {
+      Lp50 = 85.9;
+      C1 = 2.0;
+      C2 = 2.0;
+    }
+    else if (h1 > 10.0)
+    {
+      Lp50 = 76.3;
+      C1 = 2.0;
+      C2 = 1.0;
+    }
+    else if (h2 > 10.0)
+    {
+      Lp50 = 76.3;
+      C1 = 1.0;
+      C2 = 2.0;
+    }
+    else // both antenna heights below 10 metres
+    {
+      Lp50 = 66.7;
+      C1 = 1.0;
+      C2 = 1.0;
+    } // end if
+
+    Lp50 += 4.0f*_10log10f(d) + 2.0f*_10log10f(f) - C1*_10log10f(h1) - C2*_10log10f(h2);
+  /*}
+  else
+  {
+    fprintf(stderr,"Parameter error: Egli path loss model f=%6.2f h1=%6.2f h2=%6.2f d=%6.2f\n", f, h1, h2, d);
+    exit(EXIT_FAILURE);
+  }*/
+
+  return(Lp50);
+
+}
+

models/egli.hh

@@ -0,0 +1,6 @@
+#ifndef _EGLI_HH_
+#define _EGLI_HH_
+
+double EgliPathLoss(float f, float h1, float h2, float d);
+
+#endif /* _EGLI_HH_ */

models/ericsson.cc

@@ -10,12 +10,12 @@ double EricssonpathLoss(float f, float TxH, float RxH, float d, int mode)
 	// Urban 
 	double a0 = 36.2, a1 = 30.2, a2 = -12, a3 = 0.1;
 
-	if (f < 150 || f > 1900) {
+/*	if (f < 150 || f > 1900) {
 		fprintf
 		    (stderr,"Error: Ericsson9999 model frequency range 150-1900MHz\n");
 		exit(EXIT_FAILURE);
 	}
-
+*/
 	if (mode == 2) {	// Suburban / Med loss
 		a0 = 43.2;
 		a1 = 68.93;
@@ -25,8 +25,7 @@ double EricssonpathLoss(float f, float TxH, float RxH, float d, int mode)
 		a1 = 100.6;
 	}
 	double g1 = 3.2 * (log10(11.75 * RxH) * log10(11.75 * RxH));
-	double g2 = (44.49 * log10(f)) - 4.78 * ((log10(f) * log10(f)));
+	double g2 = 44.49 * log10(f) - 4.78 * (log10(f) * log10(f));
 
-	return a0 + a1 * log10(d) + a2 * log10(TxH) +
+	return a0 + a1 * log10(d) + a2 * log10(TxH) + a3 * log10(TxH) * log10(d) - g1 + g2;
-	    a3 * log10(TxH) * log10(d) - g1 + g2;
 }

models/los.cc

@@ -10,6 +10,7 @@
 #include "itwom3.0.hh"
 #include "sui.hh"
 #include "pel.hh"
+#include "egli.hh"
 #include <pthread.h>
 
 #define NUM_SECTIONS 4
@@ -304,8 +305,9 @@ void PlotPropPath(struct site source, struct site destination,
 	    xmtr_alt, dest_alt, xmtr_alt2, dest_alt2,
 	    cos_rcvr_angle, cos_test_angle = 0.0, test_alt,
 	    elevation = 0.0, distance = 0.0, four_thirds_earth,
-	    field_strength = 0.0, rxp, dBm, dkm, diffloss;
+	    field_strength = 0.0, rxp, dBm, diffloss;
 	struct site temp;
+	float dkm;
 
 	ReadPath(source, destination);
 
@@ -445,7 +447,6 @@ void PlotPropPath(struct site source, struct site destination,
 
 			dkm = (elev[1] * elev[0]) / 1000;	// km
 
-
 			switch (propmodel) {
 			case 1:
 				// Longley Rice ITM
@@ -520,12 +521,14 @@ void PlotPropPath(struct site source, struct site destination,
 						      METERS_PER_FOOT), dkm,
 						     pmenv);
 				break;
-
 			case 10:
 				// Plane earth
 				loss =	PlaneEarthLoss(dkm, source.alt * METERS_PER_FOOT, (path.elevation[y] * METERS_PER_FOOT) + (destination.alt * METERS_PER_FOOT));
 				break;
-
+			case 11:
+				// Egli VHF/UHF
+				loss = EgliPathLoss(LR.frq_mhz, source.alt * METERS_PER_FOOT, (path.elevation[y] * METERS_PER_FOOT) + (destination.alt * METERS_PER_FOOT),dkm);
+				break;
 			default:
 				point_to_point_ITM(source.alt * METERS_PER_FOOT,
 						   destination.alt *

models/sui.cc

@@ -2,49 +2,45 @@
 #include <stdlib.h>
 #include <math.h>
 
-double SUIpathLoss(float f, float TxH, float RxH, float d, int mode)
+double SUIpathLoss(double f, double TxH, double RxH, double d, int mode)
 {
 	/*
-	   f = Frequency (MHz)
+	   f = Frequency (MHz) 1900 to 11000
 	   TxH =  Transmitter height (m)
 	   RxH = Receiver height (m)
 	   d = distance (km)
-	   mode A1 = Hilly + trees
+	   mode A1 = URBAN / OBSTRUCTED
-	   mode B2 = Flat + trees OR hilly + light foliage
+	   mode B2 = SUBURBAN / PARTIALLY OBSTRUCTED
-	   mode C3 = Flat + light foliage
+	   mode C3 = RURAL / OPEN
 	   http://www.cl.cam.ac.uk/research/dtg/lce-pub/public/vsa23/VTC05_Empirical.pdf
 	 */
-	d = d * 1000;		// km to m
+	d = d * 1000.0;		// km to m
-	if (f < 1900 || f > 11000) {
+	RxH = RxH * 1000.0; // Correction factor for CPE units.
-		fprintf(stderr,"Error: SUI model frequency range 1.9-11GHz\n");
+
-		exit(EXIT_FAILURE);
+	// Urban (A1) is default
-	}
-	// Terrain mode A is default
 	double a = 4.6;
 	double b = 0.0075;
 	double c = 12.6;
-	double s = 10.6;	// Optional fading value
+	double s = 0.0; // Optional fading value. Max 10.6dB
-	int XhCF = -10.8;
+	double XhCF = -10.8;
 
-	if (mode == 2) {
+	if (mode == 2) { // Suburban
 		a = 4.0;
 		b = 0.0065;
 		c = 17.1;
-		s = 6;		// average
 	}
-	if (mode == 3) {
+	if (mode == 3) { // Rural
 		a = 3.6;
 		b = 0.005;
 		c = 20;
-		s = 3;		// Optimistic
 		XhCF = -20;
 	}
 	double d0 = 100;
-	double A = 20 * log10((4 * M_PI * d0) / (300 / f));
+	double A = 20 * log10((4 * M_PI * d0) / (300.0 / f));
 	double y = a - (b * TxH) + (c / TxH);
-	double Xf = 6 * log10(f / 2000);
+	//Correction factors
+	double Xf = 6.0 * log10(f / 2000);
 	double Xh = XhCF * log10(RxH / 2000);
 
-	//return A + (10 * y * log10(d / d0)) + Xf + Xh + s;
+	return A + (10 * y * log10(d / d0)) + Xf + Xh + s;
-	return A + (10 * y) * (log10(d / d0)) + Xf + Xh + s;
 }

models/sui.hh

@@ -1,6 +1,6 @@
 #ifndef _SUI_HH_
 #define _SUI_HH_
 
-double SUIpathLoss(float f, float TxH, float RxH, float d, int mode);
+double SUIpathLoss(double f, double TxH, double RxH, double d, int mode);
 
 #endif /* _SUI_HH_ */

models/test.cc

@@ -0,0 +1,117 @@
+#include <stdio.h>
+#include <stdlib.h>
+#include <math.h>
+#include <sys/stat.h>
+/*
+*  Propagation model test script for signal server
+*  Requires gnuplot
+*  Compile: gcc -Wall -o test test.cc sui.cc cost.cc ecc33.cc ericsson.cc fspl.cc egli.cc hata.cc -lm
+*  Test 850Mhz: ./test 850
+*
+*  This program is free software; you can redistribute it and/or modify it   *
+*  under the terms of the GNU General Public License as published by the     *
+*  Free Software Foundation; either version 2 of the License or any later    *
+*  version.                                                                  *
+*                                                                            *
+*  This program is distributed in the hope that it will useful, but WITHOUT  *
+*  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or     *
+*  FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License     *
+*  for more details.                                                         *
+*                                                                            *
+\****************************************************************************/
+
+
+extern double EgliPathLoss(float f, float TxH, float RxH, float d);
+extern double SUIpathLoss(double f, double TxH, double RxH, double d, int mode);
+extern double COST231pathLoss(float f, float TxH, float RxH, float d, int mode);
+extern double ECC33pathLoss(float f, float TxH, float RxH, float d, int mode);
+extern double EricssonpathLoss(float f, float TxH, float RxH, float d, int mode);
+extern double FSPLpathLoss(float f, float d);
+extern double HATApathLoss(float f, float TxH, float RxH, float d, int mode);
+extern void point_to_point_ITM(double tht_m, double rht_m, double eps_dielect,
+                        double sgm_conductivity, double eno_ns_surfref,
+                        double frq_mhz, int radio_climate, int pol,
+                        double conf, double rel, double &dbloss, char *strmode,
+                        int &errnum);
+
+int main(int argc, char** argv)
+{
+double a = 0;
+double f = atof(argv[1]);
+double r = 5.0;
+float TxH = 30.0;
+float RxH = 2.0;
+
+mkdir("tests", S_IRWXU | S_IRWXG | S_IROTH | S_IXOTH);
+
+FILE * fh;
+
+/*fh = fopen("tests/ITM","w");
+for(float d = 0.1; d <= r; d=d+0.2){
+    point_to_point_ITM(TxH,RxH,15.0,0.005,301.0,f,5,1,0.5,0.5,a,"",errno);
+    fprintf(fh,"%.1f\t%.1f\n",d,a);
+}
+fclose(fh);
+*/
+fh = fopen("tests/SUI.1","w");
+for(float d = 0.1; d <= r; d=d+0.1){
+    a = SUIpathLoss(f, TxH, RxH, d,1);
+    fprintf(fh,"%.1f\t%.1f\n",d,a);
+}
+fclose(fh);
+
+fh = fopen("tests/COST231.1","w");
+for(float d = 0.1; d <= r; d=d+0.1){
+    a = COST231pathLoss(f, TxH, RxH, d,1);
+    fprintf(fh,"%.1f\t%.1f\n",d,a);
+}
+fclose(fh);
+
+
+fh = fopen("tests/ECC33.1","w");
+for(float d = 0.1; d <= r; d=d+0.1){
+    a = ECC33pathLoss(f, TxH, RxH, d,1);
+    fprintf(fh,"%.1f\t%.1f\n",d,a);
+}
+fclose(fh);
+
+fh = fopen("tests/Ericsson9999.1","w");
+for(float d = 0.1; d <= r; d=d+0.1){
+    a = EricssonpathLoss(f, TxH, RxH, d,1);
+    fprintf(fh,"%.1f\t%.1f\n",d,a);
+}
+fclose(fh);
+
+fh = fopen("tests/FSPL","w");
+for(float d = 0.1; d <= r; d=d+0.1){
+    a = FSPLpathLoss(f, d);
+    fprintf(fh,"%.1f\t%.1f\n",d,a);
+}
+fclose(fh);
+
+fh = fopen("tests/Hata.1","w");
+for(float d = 0.1; d <= r; d=d+0.1){
+    a = HATApathLoss(f, TxH, RxH, d, 1);
+    fprintf(fh,"%.1f\t%.1f\n",d,a);
+}
+fclose(fh);
+
+fh = fopen("tests/Egli.VHF-UHF","w");
+for(float d = 0.1; d <= r; d=d+0.1){
+    a = EgliPathLoss(f, TxH, RxH, d);
+    fprintf(fh,"%.1f\t%.1f\n",d,a);
+}
+fclose(fh);
+
+
+fh = fopen("tests/gnuplot.plt","w");
+fprintf(fh,"set terminal jpeg size 800,600\nset output '%.0fMHz_propagation_models.jpg'\nset title '%.0fMHz path loss by propagation model - CloudRF.com'\n",f,f);
+fprintf(fh,"set key right bottom box\nset style line 1\nset grid\nset xlabel 'Distance KM'\nset ylabel 'Path Loss dB'\n");
+fprintf(fh,"plot 'FSPL' with lines, 'Hata.1' with lines, 'COST231.1' with lines,'ECC33.1' with lines,'Ericsson9999.1' with lines,'SUI.1' with lines,'Egli.VHF-UHF' with lines");
+fclose(fh);
+
+system("cd tests && gnuplot gnuplot.plt");
+
+return(0);
+}
+

outputs.cc

@@ -1534,13 +1534,15 @@ void PathReport(struct site source, struct site destination, char *name,
 				free_space_loss);
 		}
 
-		if((loss*1.1) < free_space_loss){
-			fprintf(stderr,"Model error! Computed loss of %.1fdB is greater than free space loss of %.1fdB. Check your inuts for model %d\n",loss,free_space_loss,propmodel);
-			return;
-		}
-
 		fprintf(fd2, "Computed path loss: %.2f dB\n", loss);
 
+
+                if((loss*1.5) < free_space_loss){
+			fprintf(fd2,"Model error! Computed loss of %.1fdB is greater than free space loss of %.1fdB. Check your inuts for model %d\n",loss,free_space_loss,propmodel);
+                        fprintf(stderr,"Model error! Computed loss of %.1fdB is greater than free space loss of %.1fdB. Check your inuts for model %d\n",loss,free_space_loss,propmodel);
+                        return;
+                }
+
 		if (free_space_loss != 0.0)
 			fprintf(fd2,
 				"Attenuation due to terrain shielding: %.2f dB\n",