From aca19119fd7790b083e6777d7ba29f8da0f83c06 Mon Sep 17 00:00:00 2001 From: Harald Welte Date: Thu, 20 Aug 2015 22:32:35 +0200 Subject: [PATCH 1/8] Add support for transmission via UHD devices this is currently only tested with an USRP2, but should work similarly with other UHD devices. --- README.md | 13 +++++++- gps-sdr-sim-uhd.py | 74 ++++++++++++++++++++++++++++++++++++++++++++++ 2 files changed, 86 insertions(+), 1 deletion(-) create mode 100755 gps-sdr-sim-uhd.py diff --git a/README.md b/README.md index 26d93e9..60d4649 100644 --- a/README.md +++ b/README.md @@ -38,7 +38,11 @@ Doppler for the GPS satellites in view. This simulated range data is then used to generate the digitized I/Q samples for the GPS signal. The bladeRF command line interface requires I/Q pairs stored as signed -16-bit integers, while the hackrf_transfer supports signed bytes. +16-bit integers, while the hackrf_transfer and gps-sdr-sim-uhd.py +supports signed bytes. + +HackRF + bladeRF require 2.6 MHz sample rate, while the USRP2 requires +2.5 MHz (an even integral decimator of 100 MHz). ``` Usage: gps-sdr-sim [options] @@ -92,6 +96,13 @@ For the HackRF: > hackrf_transfer -t gpssim.bin -f 1575420000 -s 2600000 -a 1 -x 0 ``` +For UHD supported devices (tested with USRP2 only): + +``` +> gps-sdr-sim-uhd.py -t gpssim.bin -s 2500000 -x 0 +``` + + ### License Copyright © 2015 Takuji Ebinuma diff --git a/gps-sdr-sim-uhd.py b/gps-sdr-sim-uhd.py new file mode 100755 index 0000000..06da2ce --- /dev/null +++ b/gps-sdr-sim-uhd.py @@ -0,0 +1,74 @@ +#!/usr/bin/env python +# a small script to transmit simulated GPS samples via UHD +# (C) 2015 by Harald Welte +# Licensed under the MIT License (see LICENSE) + +from gnuradio import blocks +from gnuradio import eng_notation +from gnuradio import gr +from gnuradio import uhd +from gnuradio.eng_option import eng_option +from gnuradio.filter import firdes +from optparse import OptionParser +import time + +class top_block(gr.top_block): + + def __init__(self, options): + gr.top_block.__init__(self, "GPS-SDR-SIM") + + ################################################## + # Blocks + ################################################## + self.uhd_usrp_sink = uhd.usrp_sink( + ",".join(("", "")), + uhd.stream_args( + cpu_format="fc32", + channels=range(1), + ), + ) + self.uhd_usrp_sink.set_samp_rate(options.sample_rate) + self.uhd_usrp_sink.set_center_freq(options.frequency, 0) + self.uhd_usrp_sink.set_gain(options.gain, 0) + + # a file source for the file generated by the gps-sdr-sim + self.blocks_file_source = blocks.file_source(gr.sizeof_char*1, options.filename, True) + + # convert from signed bytes to short + self.blocks_char_to_short = blocks.char_to_short(1) + + # convert from interleaved short to complex values + self.blocks_interleaved_short_to_complex = blocks.interleaved_short_to_complex(False, False) + + # establish the connections + self.connect((self.blocks_file_source, 0), (self.blocks_char_to_short, 0)) + self.connect((self.blocks_char_to_short, 0), (self.blocks_interleaved_short_to_complex, 0)) + self.connect((self.blocks_interleaved_short_to_complex, 0), (self.uhd_usrp_sink, 0)) + +def get_options(): + parser = OptionParser(option_class=eng_option) + parser.add_option("-x", "--gain", type="eng_float", default=0, + help="set transmitter gain [default=0]") + parser.add_option("-f", "--frequency", type="eng_float", default=1575420000, + help="set transmit frequency [default=1575420000]") + # On USRP2, the sample rate should lead to an even decimator + # based on the 100 MHz clock. At 2.5 MHz, we end up with 40 + parser.add_option("-s", "--sample-rate", type="eng_float", default=2500000, + help="set sample rate [default=2500000]") + parser.add_option("-t", "--filename", type="string", default="gpssim.bin", + help="set output file name [default=gpssim.bin]") + + (options, args) = parser.parse_args() + if len(args) != 0: + parser.print_help() + raise SystemExit, 1 + + return (options) + +if __name__ == '__main__': + (options) = get_options() + tb = top_block(options) + tb.start() + raw_input('Press Enter to quit: ') + tb.stop() + tb.wait() From 22f01222750cfed59b38b3c61200c3ad1c468cd4 Mon Sep 17 00:00:00 2001 From: Harald Welte Date: Thu, 20 Aug 2015 22:30:54 +0200 Subject: [PATCH 2/8] gpssim.c: Use 'const' keyword for read-only function input arguments It is generally good coding practise to use the const qualifier on arguments that are pointers to input-only arguments of functions. --- gpssim.c | 24 ++++++++++++------------ 1 file changed, 12 insertions(+), 12 deletions(-) diff --git a/gpssim.c b/gpssim.c index 85bf505..5c40e57 100644 --- a/gpssim.c +++ b/gpssim.c @@ -247,7 +247,7 @@ void codegen(int *ca, int prn) return; } -void date2gps(datetime_t *t, gpstime_t *g) +void date2gps(const datetime_t *t, gpstime_t *g) { int doy[12] = {0,31,59,90,120,151,181,212,243,273,304,334}; int ye; @@ -272,7 +272,7 @@ void date2gps(datetime_t *t, gpstime_t *g) return; } -void xyz2llh(double *xyz, double *llh) +void xyz2llh(const double *xyz, double *llh) { double a,eps,e,e2; double x,y,z; @@ -312,7 +312,7 @@ void xyz2llh(double *xyz, double *llh) return; } -void llh2xyz(double *llh, double *xyz) +void llh2xyz(const const double *llh, double *xyz) { double n; double a; @@ -346,7 +346,7 @@ void llh2xyz(double *llh, double *xyz) return; } -void ltcmat(double *llh, double t[3][3]) +void ltcmat(const double *llh, double t[3][3]) { double slat, clat; double slon, clon; @@ -369,7 +369,7 @@ void ltcmat(double *llh, double t[3][3]) return; } -void ecef2neu(double *xyz, double t[3][3], double *neu) +void ecef2neu(const double *xyz, double t[3][3], double *neu) { neu[0] = t[0][0]*xyz[0] + t[0][1]*xyz[1] + t[0][2]*xyz[2]; neu[1] = t[1][0]*xyz[0] + t[1][1]*xyz[1] + t[1][2]*xyz[2]; @@ -378,7 +378,7 @@ void ecef2neu(double *xyz, double t[3][3], double *neu) return; } -void neu2azel(double *azel, double *neu) +void neu2azel(double *azel, const double *neu) { double ne; @@ -498,7 +498,7 @@ void satpos(ephem_t eph, gpstime_t g, double *pos, double *vel, double *clk) return; } -void eph2sbf(ephem_t eph, unsigned long sbf[5][10]) +void eph2sbf(const ephem_t eph, unsigned long sbf[5][10]) { unsigned long wn; unsigned long toe; @@ -722,7 +722,7 @@ int checkExpDesignator(char *str, int len) return(n); } -int readRinexNav(ephem_t eph[], char *fname) +int readRinexNav(ephem_t eph[], const char *fname) { FILE *fp; int nsat; @@ -963,7 +963,7 @@ int readRinexNav(ephem_t eph[], char *fname) return(nsat); } -void subVect(double *y, double *x1, double *x2) +void subVect(double *y, const double *x1, const double *x2) { y[0] = x1[0]-x2[0]; y[1] = x1[1]-x2[1]; @@ -972,12 +972,12 @@ void subVect(double *y, double *x1, double *x2) return; } -double normVect(double *x) +double normVect(const double *x) { return(sqrt(x[0]*x[0]+x[1]*x[1]+x[2]*x[2])); } -double dotProd(double *x1, double *x2) +double dotProd(const double *x1, const double *x2) { return(x1[0]*x2[0]+x1[1]*x2[1]+x1[2]*x2[2]); } @@ -1060,7 +1060,7 @@ void computeCodePhase(channel_t *chan, range_t rho0, range_t rho1, double dt) return; } -int readUserMotion(double xyz[USER_MOTION_SIZE][3], char *filename) +int readUserMotion(double xyz[USER_MOTION_SIZE][3], const char *filename) { FILE *fp; int numd; From 78dd5d7266e1d787b68f35efc3eabf878fb76771 Mon Sep 17 00:00:00 2001 From: Harald Welte Date: Thu, 20 Aug 2015 22:31:50 +0200 Subject: [PATCH 3/8] Makefile: Add -Wall to generate more compiler warnings --- Makefile | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Makefile b/Makefile index 96aa71b..b741004 100644 --- a/Makefile +++ b/Makefile @@ -5,7 +5,7 @@ all: gps-sdr-sim SHELL=/bin/bash CC=gcc -CFLAGS=-fopenmp -O3 +CFLAGS=-fopenmp -O3 -Wall LDFLAGS=-lm -fopenmp gps-sdr-sim: gpssim.o From 00202f91dd05a9d8ea84d218593212367f9a01b0 Mon Sep 17 00:00:00 2001 From: Harald Welte Date: Fri, 21 Aug 2015 00:00:06 +0200 Subject: [PATCH 4/8] gpssim: Add docbook-style documentation to structures + functions I felt that the current code could deal with a lot more comments, particularly with somebody not familiar with all the very-few-letter mnemonics of th pseudorange equation. It is my hope that by adding these comments, somebody with a basic understanding of GPS can read through the code more easily. --- gpssim.c | 203 ++++++++++++++++++++++++++++++++++++++++++------------- 1 file changed, 156 insertions(+), 47 deletions(-) diff --git a/gpssim.c b/gpssim.c index 5c40e57..2bc4bed 100644 --- a/gpssim.c +++ b/gpssim.c @@ -18,14 +18,21 @@ typedef int bool; #define false 0 #endif +/*! \brief Maximum length of a line in a text file (RINEX, motion) */ #define MAX_CHAR (100) +/*! \brief Maximum number of satellites in RINEX file */ #define MAX_SAT (32) +/*! \brief Maximum number of channels we simulate */ #define MAX_CHAN (16) +/*! \brief Maximum number of user motion waypoints */ #define USER_MOTION_SIZE (3000) // max 300 sec at 10Hz +/*! \brief Number of subframes */ #define N_SBF (51) // 6 seconds per subframe, 6 sec * 51 = 306 sec (max) + +/*! \brief Number of words */ #define N_DWRD (N_SBF*10) // 10 word per subframe #define SECONDS_IN_WEEK 604800.0 @@ -52,7 +59,9 @@ typedef int bool; #define SPEED_OF_LIGHT 2.99792458e8 #define LAMBDA_L1 0.190293672798365 +/*! \brief GPS L1 Carrier frequency */ #define CARR_FREQ (1575.42e6) +/*! \brief C/A code frequency */ #define CODE_FREQ (1.023e6) #define CARR_TO_CODE (1.0/1540.0) @@ -135,54 +144,61 @@ int cosTable512[] = { }; #endif +/*! \file gpssim.c + * \brief GPS Satellite Simulator + */ + +/*! \brief Structure representing GPS time */ typedef struct { - int week; - double sec; + int week; /*!< GPS week number (since January 1980) */ + double sec; /*!< second inside the GPS \a week */ } gpstime_t; +/*! \brief Structure repreenting UTC time */ typedef struct { - int y; - int m; - int d; - int hh; - int mm; - double sec; + int y; /*!< Calendar year */ + int m; /*!< Calendar month */ + int d; /*!< Calendar day */ + int hh; /*!< Calendar hour */ + int mm; /*!< Calendar minutes */ + double sec; /*!< Calendar seconds */ } datetime_t; +/*! \brief Structure representing ephemeris of a single satellite */ typedef struct { - int vflg; - gpstime_t toc; - gpstime_t toe; - int iodc; - int iode; - double deltan; - double cuc; - double cus; - double cic; - double cis; - double crc; - double crs; - double ecc; - double sqrta; - double m0; - double omg0; - double inc0; + int vflg; /*!< Valid Flag */ + gpstime_t toc; /*!< Time of Clock */ + gpstime_t toe; /*!< Time of Ephemeris */ + int iodc; /*!< Issue of Data, Clock */ + int iode; /*!< Isuse of Data, Ephemeris */ + double deltan; /*!< Delta-N (radians/sec) */ + double cuc; /*!< Cuc (radians) */ + double cus; /*!< Cus (radians) */ + double cic; /*!< Correction to inclination cos (radians) */ + double cis; /*!< Correction to inclination sin (radians) */ + double crc; /*!< Correction to radius cos (meters) */ + double crs; /*!< Correction to radius sin (meters) */ + double ecc; /*!< e Eccentricity */ + double sqrta; /*!< sqrt(A) (sqrt(m)) */ + double m0; /*!< Mean anamoly (radians) */ + double omg0; /*!< Longitude of the ascending node (radians) */ + double inc0; /*!< Inclination (radians) */ double aop; - double omgdot; - double idot; - double af0; - double af1; - double af2; - double tgd; + double omgdot; /*!< Omega dot (radians/s) */ + double idot; /*!< IDOT (radians/s) */ + double af0; /*!< Clock offset (seconds) */ + double af1; /*!< rate (sec/sec) */ + double af2; /*!< acceleration (sec/sec^2) */ + double tgd; /*!< Group delay L2 bias */ // Working variables follow - double n; // Mean motion (Average angular velocity) - double sq1e2; // sqrt(1-e^2) - double A; // Semi-major axis - double omgkdot; // OmegaDot-OmegaEdot + double n; /*!< Mean motion (Average angular velocity) */ + double sq1e2; /*!< sqrt(1-e^2) */ + double A; /*!< Semi-major axis */ + double omgkdot; /*!< OmegaDot-OmegaEdot */ } ephem_t; typedef struct @@ -192,22 +208,29 @@ typedef struct double rate; } range_t; +/*! \brief Structure representing a Channel */ typedef struct { - int prn; - int ca[1023]; - double f_carr,f_code; - double carr_phase,code_phase; - gpstime_t g0; - unsigned long dwrd[N_DWRD]; - int iword; - int ibit; - int icode; - int dataBit; - int codeCA; - short *iq_buff; + int prn; /*< PRN Number */ + int ca[1023]; /*< C/A Sequence */ + double f_carr; /*< Carrier frequency */ + double f_code; /*< Code frequency */ + double carr_phase; /*< Carrier phase */ + double code_phase; /*< Code phase */ + gpstime_t g0; /*!< GPS time at start */ + unsigned long dwrd[N_DWRD]; /*!< Data words of sub-frame */ + int iword; /*!< initial word */ + int ibit; /*!< initial bit */ + int icode; /*!< initial code */ + int dataBit; /*!< current data bit */ + int codeCA; /*!< current C/A code */ + short *iq_buff; /*< buffer of I/Q samples */ } channel_t; +/* !\brief generate the C/A code sequence for a given Satellite Vehicle PRN + * \param[in] prn PRN nuber of the Satellite Vehicle + * \param[out] ca Caller-allocated integer array of 1023 bytes + */ void codegen(int *ca, int prn) { int delay[] = { @@ -247,6 +270,10 @@ void codegen(int *ca, int prn) return; } +/*! \brief Convert a UTC date into a GPS date + * \param[in] t input date in UTC form + * \param[out] g output date in GPS form + */ void date2gps(const datetime_t *t, gpstime_t *g) { int doy[12] = {0,31,59,90,120,151,181,212,243,273,304,334}; @@ -272,6 +299,10 @@ void date2gps(const datetime_t *t, gpstime_t *g) return; } +/*! \brief Convert Earth-centered Earth-fixed (ECEF) into Lat/Long/Heighth + * \param[in] xyz Input Array of X, Y and Z ECEF coordinates + * \param[out] llh Output Array of Latitude, Longitude and Height + */ void xyz2llh(const double *xyz, double *llh) { double a,eps,e,e2; @@ -312,6 +343,10 @@ void xyz2llh(const double *xyz, double *llh) return; } +/*! \brief Convert Lat/Long/Height into Earth-centered Earth-fixed (ECEF) + * \param[in] llh Input Array of Latitude, Longitude and Height + * \param[out] xyz Output Array of X, Y and Z ECEF coordinates + */ void llh2xyz(const const double *llh, double *xyz) { double n; @@ -346,6 +381,10 @@ void llh2xyz(const const double *llh, double *xyz) return; } +/*! \brief Compute the intermediate matrix for LLH to ECEF + * \param[in] llh Input position in Latitude-Longitude-Height format + * \param[out] t Three-by-Three output matrix + */ void ltcmat(const double *llh, double t[3][3]) { double slat, clat; @@ -369,6 +408,11 @@ void ltcmat(const double *llh, double t[3][3]) return; } +/*! \brief Convert Earth-centered Earth-Fixed to ? + * \param[in] xyz Input position as vector in ECEF format + * \param[in] t Intermediate matrix computed by \ref ltcmat + * \param[out] neu Output position as North-East-Up format + */ void ecef2neu(const double *xyz, double t[3][3], double *neu) { neu[0] = t[0][0]*xyz[0] + t[0][1]*xyz[1] + t[0][2]*xyz[2]; @@ -378,6 +422,10 @@ void ecef2neu(const double *xyz, double t[3][3], double *neu) return; } +/*! \brief Convert North-Eeast-Up to Azimuth + Elevation + * \param[in] neu Input position in North-East-Up format + * \param[out] azel Output array of azimuth + elevation as double + */ void neu2azel(double *azel, const double *neu) { double ne; @@ -392,6 +440,13 @@ void neu2azel(double *azel, const double *neu) return; } +/*! \brief Compute Satellite position, velocity and clock at given time + * \param[in] eph Ephemeris data of the satellite + * \param[in] g GPS time at which position is to be computed + * \param[out] pos Computed position (vector) + * \param[out] vel Computed velociy (vector) + * \param[clk] clk Computed clock + */ void satpos(ephem_t eph, gpstime_t g, double *pos, double *vel, double *clk) { // Computing Satellite Velocity using the Broadcast Ephemeris @@ -498,6 +553,10 @@ void satpos(ephem_t eph, gpstime_t g, double *pos, double *vel, double *clk) return; } +/*! \brief Compute Subframe from Ephemeris + * \param[in] eph Ephemeris of given SV + * \param[out] sbf Array of five sub-frames, 10 long words each + */ void eph2sbf(const ephem_t eph, unsigned long sbf[5][10]) { unsigned long wn; @@ -624,6 +683,10 @@ void eph2sbf(const ephem_t eph, unsigned long sbf[5][10]) return; } +/*! \brief Count number of bits set to 1 + * \param[in] v long word in whihc bits are counted + * \returns Count of bits set to 1 + */ unsigned long countBits(unsigned long v) { unsigned long c; @@ -641,6 +704,11 @@ unsigned long countBits(unsigned long v) return(c); } +/*! \brief Compute the Checksum for one given word of a subframe + * \param[in] source The input data + * \param[in] nib Does this word contain non-information-bearing bits? + * \returns Computed Checksum + */ unsigned long computeChecksum(unsigned long source, int nib) { /* @@ -706,6 +774,11 @@ unsigned long computeChecksum(unsigned long source, int nib) return(D); } +/*! \brief Replace all 'E' exponential designators to 'D' + * \param str String in which all occurrences of 'E' are replaced with * 'D' + * \param len Length of input string in bytes + * \returns Number of characters replaced + */ int checkExpDesignator(char *str, int len) { int i,n=0; @@ -722,6 +795,11 @@ int checkExpDesignator(char *str, int len) return(n); } +/*! \brief Read Ephemersi data from the RINEX Navigation file */ +/* \param[out] eph Array of Output SV ephemeris data + * \param[in] fname File name of the RINEX file + * \returns Number of SV found in the file, -1 on error + */ int readRinexNav(ephem_t eph[], const char *fname) { FILE *fp; @@ -963,6 +1041,11 @@ int readRinexNav(ephem_t eph[], const char *fname) return(nsat); } +/*! \brief Subtract two vectors of double + * \param[out] y Result of subtraction + * \param[in] x1 Minuend of subtracion + * \param[in] x2 Subtrahend of subtracion + */ void subVect(double *y, const double *x1, const double *x2) { y[0] = x1[0]-x2[0]; @@ -972,16 +1055,31 @@ void subVect(double *y, const double *x1, const double *x2) return; } +/*! \brief Compute Norm of Vector + * \param[in] x Input vector + * \returns Length (Norm) of the input vector + */ double normVect(const double *x) { return(sqrt(x[0]*x[0]+x[1]*x[1]+x[2]*x[2])); } +/*! \brief Compute dot-product of two vectors + * \param[in] x1 First multiplicand + * \param[in] x2 Second multiplicand + * \returns Dot-product of both multiplicands + */ double dotProd(const double *x1, const double *x2) { return(x1[0]*x2[0]+x1[1]*x2[1]+x1[2]*x2[2]); } +/*! \brief Compute range between a satellite and the receiver + * \param[out] rho The computed range + * \param[in] eph Ephemeris data of the satellite + * \param[in] g GPS time at time of receiving the signal + * \param[in] xyz position of the receiver + */ void computeRange(range_t *rho, ephem_t eph, gpstime_t g, double xyz[]) { double pos[3],vel[3],clk[2]; @@ -1027,6 +1125,12 @@ void computeRange(range_t *rho, ephem_t eph, gpstime_t g, double xyz[]) return; } +/*! \brief Compute the code phase for a given channel (satellite) + * \param chan Channel on which we operate (is updated) + * \param[in] rho0 Range at start of interval + * \param[in] rho1 Current range, after \a dt has expired + * \param[in dt delta-t (time difference) in seconds + */ void computeCodePhase(channel_t *chan, range_t rho0, range_t rho1, double dt) { double ms; @@ -1060,6 +1164,11 @@ void computeCodePhase(channel_t *chan, range_t rho0, range_t rho1, double dt) return; } +/*! \brief Read the list of user motions from the input file + * \param[out] xyz Output array of ECEF vectors for user motion + * \param[[in] filename File name of the text input file + * \returns Number of user data motion records read, -1 on error + */ int readUserMotion(double xyz[USER_MOTION_SIZE][3], const char *filename) { FILE *fp; From f1520f0cce1991be09ad3b9edf732ab86a4161fa Mon Sep 17 00:00:00 2001 From: Harald Welte Date: Fri, 21 Aug 2015 00:10:52 +0200 Subject: [PATCH 5/8] Use #defines instead of WGS84 magic numbers occurring repeatedly --- gpssim.c | 11 +++++++---- 1 file changed, 7 insertions(+), 4 deletions(-) diff --git a/gpssim.c b/gpssim.c index 2bc4bed..a1a6d76 100644 --- a/gpssim.c +++ b/gpssim.c @@ -54,6 +54,9 @@ typedef int bool; #define OMEGA_EARTH 7.2921151467e-5 #define PI 3.1415926535898 +#define WGS84_RADIUS 6378137.0 +#define WGS84_ECCENTRICITY 0.0818191908426 + #define R2D 57.2957795131 #define SPEED_OF_LIGHT 2.99792458e8 @@ -309,8 +312,8 @@ void xyz2llh(const double *xyz, double *llh) double x,y,z; double rho2,dz,zdz,nh,slat,n,dz_new; - a = 6378137.0; - e = 0.0818191908426; + a = WGS84_RADIUS; + e = WGS84_ECCENTRICITY; eps = 1.0e-3; e2 = e*e; @@ -360,8 +363,8 @@ void llh2xyz(const const double *llh, double *xyz) double d,nph; double tmp; - a = 6378137.0; - e = 0.0818191908426; + a = WGS84_RADIUS; + e = WGS84_ECCENTRICITY; e2 = e*e; clat = cos(llh[0]); From 4282f2b63e0464d4311f9d9b4363bf00556dc8bf Mon Sep 17 00:00:00 2001 From: Harald Welte Date: Fri, 21 Aug 2015 00:59:23 +0200 Subject: [PATCH 6/8] rename checkExpDesignator to replaceExpDesignator 'check' is somewhat of a misnomer, as the function doens't just check something, but it actively replaces something. --- gpssim.c | 48 ++++++++++++++++++++++++------------------------ 1 file changed, 24 insertions(+), 24 deletions(-) diff --git a/gpssim.c b/gpssim.c index a1a6d76..e60a465 100644 --- a/gpssim.c +++ b/gpssim.c @@ -782,7 +782,7 @@ unsigned long computeChecksum(unsigned long source, int nib) * \param len Length of input string in bytes * \returns Number of characters replaced */ -int checkExpDesignator(char *str, int len) +int replaceExpDesignator(char *str, int len) { int i,n=0; @@ -883,17 +883,17 @@ int readRinexNav(ephem_t eph[], const char *fname) strncpy(tmp, str+22, 19); tmp[19] = 0; - checkExpDesignator(tmp, 19); // tmp[15]='E'; + replaceExpDesignator(tmp, 19); // tmp[15]='E'; eph[sv].af0 = atof(tmp); strncpy(tmp, str+41, 19); tmp[19] = 0; - checkExpDesignator(tmp, 19); + replaceExpDesignator(tmp, 19); eph[sv].af1 = atof(tmp); strncpy(tmp, str+60, 19); tmp[19] = 0; - checkExpDesignator(tmp, 19); + replaceExpDesignator(tmp, 19); eph[sv].af2 = atof(tmp); // BROADCAST ORBIT - 1 @@ -902,22 +902,22 @@ int readRinexNav(ephem_t eph[], const char *fname) strncpy(tmp, str+3, 19); tmp[19] = 0; - checkExpDesignator(tmp, 19); + replaceExpDesignator(tmp, 19); eph[sv].iode = (int)atof(tmp); strncpy(tmp, str+22, 19); tmp[19] = 0; - checkExpDesignator(tmp, 19); + replaceExpDesignator(tmp, 19); eph[sv].crs = atof(tmp); strncpy(tmp, str+41, 19); tmp[19] = 0; - checkExpDesignator(tmp, 19); + replaceExpDesignator(tmp, 19); eph[sv].deltan = atof(tmp); strncpy(tmp, str+60, 19); tmp[19] = 0; - checkExpDesignator(tmp, 19); + replaceExpDesignator(tmp, 19); eph[sv].m0 = atof(tmp); // BROADCAST ORBIT - 2 @@ -926,22 +926,22 @@ int readRinexNav(ephem_t eph[], const char *fname) strncpy(tmp, str+3, 19); tmp[19] = 0; - checkExpDesignator(tmp, 19); + replaceExpDesignator(tmp, 19); eph[sv].cuc = atof(tmp); strncpy(tmp, str+22, 19); tmp[19] = 0; - checkExpDesignator(tmp, 19); + replaceExpDesignator(tmp, 19); eph[sv].ecc = atof(tmp); strncpy(tmp, str+41, 19); tmp[19] = 0; - checkExpDesignator(tmp, 19); + replaceExpDesignator(tmp, 19); eph[sv].cus = atof(tmp); strncpy(tmp, str+60, 19); tmp[19] = 0; - checkExpDesignator(tmp, 19); + replaceExpDesignator(tmp, 19); eph[sv].sqrta = atof(tmp); // BROADCAST ORBIT - 3 @@ -950,22 +950,22 @@ int readRinexNav(ephem_t eph[], const char *fname) strncpy(tmp, str+3, 19); tmp[19] = 0; - checkExpDesignator(tmp, 19); + replaceExpDesignator(tmp, 19); eph[sv].toe.sec = atof(tmp); strncpy(tmp, str+22, 19); tmp[19] = 0; - checkExpDesignator(tmp, 19); + replaceExpDesignator(tmp, 19); eph[sv].cic = atof(tmp); strncpy(tmp, str+41, 19); tmp[19] = 0; - checkExpDesignator(tmp, 19); + replaceExpDesignator(tmp, 19); eph[sv].omg0 = atof(tmp); strncpy(tmp, str+60, 19); tmp[19] = 0; - checkExpDesignator(tmp, 19); + replaceExpDesignator(tmp, 19); eph[sv].cis = atof(tmp); // BROADCAST ORBIT - 4 @@ -974,22 +974,22 @@ int readRinexNav(ephem_t eph[], const char *fname) strncpy(tmp, str+3, 19); tmp[19] = 0; - checkExpDesignator(tmp, 19); + replaceExpDesignator(tmp, 19); eph[sv].inc0 = atof(tmp); strncpy(tmp, str+22, 19); tmp[19] = 0; - checkExpDesignator(tmp, 19); + replaceExpDesignator(tmp, 19); eph[sv].crc = atof(tmp); strncpy(tmp, str+41, 19); tmp[19] = 0; - checkExpDesignator(tmp, 19); + replaceExpDesignator(tmp, 19); eph[sv].aop = atof(tmp); strncpy(tmp, str+60, 19); tmp[19] = 0; - checkExpDesignator(tmp, 19); + replaceExpDesignator(tmp, 19); eph[sv].omgdot = atof(tmp); // BROADCAST ORBIT - 5 @@ -998,12 +998,12 @@ int readRinexNav(ephem_t eph[], const char *fname) strncpy(tmp, str+3, 19); tmp[19] = 0; - checkExpDesignator(tmp, 19); + replaceExpDesignator(tmp, 19); eph[sv].idot = atof(tmp); strncpy(tmp, str+41, 19); tmp[19] = 0; - checkExpDesignator(tmp, 19); + replaceExpDesignator(tmp, 19); eph[sv].toe.week = (int)atof(tmp); // BROADCAST ORBIT - 6 @@ -1012,12 +1012,12 @@ int readRinexNav(ephem_t eph[], const char *fname) strncpy(tmp, str+41, 19); tmp[19] = 0; - checkExpDesignator(tmp, 19); + replaceExpDesignator(tmp, 19); eph[sv].tgd = atof(tmp); strncpy(tmp, str+60, 19); tmp[19] = 0; - checkExpDesignator(tmp, 19); + replaceExpDesignator(tmp, 19); eph[sv].iodc = (int)atof(tmp); // BROADCAST ORBIT - 7 From 9081cad620b1b6af908205252bc8e426d4660378 Mon Sep 17 00:00:00 2001 From: Harald Welte Date: Fri, 21 Aug 2015 01:10:53 +0200 Subject: [PATCH 7/8] Introduce N_DWRD_SBF as symbolic constant for 10 dwords-per-subframe --- gpssim.c | 19 +++++++++++-------- 1 file changed, 11 insertions(+), 8 deletions(-) diff --git a/gpssim.c b/gpssim.c index e60a465..3815d0a 100644 --- a/gpssim.c +++ b/gpssim.c @@ -32,8 +32,11 @@ typedef int bool; /*! \brief Number of subframes */ #define N_SBF (51) // 6 seconds per subframe, 6 sec * 51 = 306 sec (max) +/*! \brief Number of words per subframe */ +#define N_DWRD_SBF (10) // 10 word per subframe + /*! \brief Number of words */ -#define N_DWRD (N_SBF*10) // 10 word per subframe +#define N_DWRD (N_SBF*N_DWRD_SBF) // 10 word per subframe #define SECONDS_IN_WEEK 604800.0 #define SECONDS_IN_HALF_WEEK 302400.0 @@ -242,13 +245,13 @@ void codegen(int *ca, int prn) 473, 474, 509, 512, 513, 514, 515, 516, 859, 860, 861, 862}; - int g1[1023],g2[1023],r1[10],r2[10],c1,c2; + int g1[1023], g2[1023], r1[N_DWRD_SBF], r2[N_DWRD_SBF], c1, c2; int i,j; if (prn<1 || prn>32) return; - for (i=0;i<10;i++) + for (i=0; i Date: Fri, 21 Aug 2015 01:14:38 +0200 Subject: [PATCH 8/8] Introduce CA_SEQ_LEN as symbolic constant for 1023 We don't need to distinguish int and double, as the int constant will be converted to double before any comparison anyway. --- gpssim.c | 21 +++++++++++++-------- 1 file changed, 13 insertions(+), 8 deletions(-) diff --git a/gpssim.c b/gpssim.c index 3815d0a..17638ce 100644 --- a/gpssim.c +++ b/gpssim.c @@ -38,6 +38,9 @@ typedef int bool; /*! \brief Number of words */ #define N_DWRD (N_SBF*N_DWRD_SBF) // 10 word per subframe +/*! \brief C/A code sequence length */ +#define CA_SEQ_LEN (1023) + #define SECONDS_IN_WEEK 604800.0 #define SECONDS_IN_HALF_WEEK 302400.0 #define SECONDS_IN_DAY 86400.0 @@ -218,7 +221,7 @@ typedef struct typedef struct { int prn; /*< PRN Number */ - int ca[1023]; /*< C/A Sequence */ + int ca[CA_SEQ_LEN]; /*< C/A Sequence */ double f_carr; /*< Carrier frequency */ double f_code; /*< Code frequency */ double carr_phase; /*< Carrier phase */ @@ -245,7 +248,9 @@ void codegen(int *ca, int prn) 473, 474, 509, 512, 513, 514, 515, 516, 859, 860, 861, 862}; - int g1[1023], g2[1023], r1[N_DWRD_SBF], r2[N_DWRD_SBF], c1, c2; + int g1[CA_SEQ_LEN], g2[CA_SEQ_LEN]; + int r1[N_DWRD_SBF], r2[N_DWRD_SBF]; + int c1, c2; int i,j; if (prn<1 || prn>32) @@ -254,7 +259,7 @@ void codegen(int *ca, int prn) for (i=0; ig0.sec)+6.0) - rho0.range/SPEED_OF_LIGHT)*1000.0; ims = (int)ms; - chan->code_phase = (ms-(double)ims)*1023.0; // in chip + chan->code_phase = (ms-(double)ims)*CA_SEQ_LEN; // in chip chan->iword = ims/600; // 1 word = 30 bits = 600 ms ims -= chan->iword*600; @@ -1633,9 +1638,9 @@ int main(int argc, char *argv[]) // Update code phase chan[i].code_phase += chan[i].f_code * delt; - if (chan[i].code_phase>=1023.0) + if (chan[i].code_phase>=CA_SEQ_LEN) { - chan[i].code_phase -= 1023.0; + chan[i].code_phase -= CA_SEQ_LEN; chan[i].icode++;