From afa374367da135226a3995b7ede987456884e6c8 Mon Sep 17 00:00:00 2001
From: joanvallve <jvallve@iri.upc.edu>
Date: Fri, 8 May 2020 18:56:52 +0200
Subject: [PATCH] debugging and not taking sbas correction error as ephemeris
 error

---
 src/ephemeris.c |  2 +-
 src/pntpos.c    | 29 ++++++++++++++++++++++-------
 2 files changed, 23 insertions(+), 8 deletions(-)

diff --git a/src/ephemeris.c b/src/ephemeris.c
index 8f62ba2..75a6d2a 100644
--- a/src/ephemeris.c
+++ b/src/ephemeris.c
@@ -563,7 +563,7 @@ static int satpos_sbas(gtime_t time, gtime_t teph, int sat, const nav_t *nav,
     if (i>=nav->sbssat.nsat) {
         trace(2,"no sbas correction for orbit: %s sat=%2d\n",time_str(time,0),sat);
         ephpos(time,teph,sat,nav,-1,rs,dts,var,svh);
-        *svh=-1;
+        *svh=-2;// JV: error code for detecting available ephemeris without sbas correction
         return 0;
     }
     /* satellite postion and clock by broadcast ephemeris */
diff --git a/src/pntpos.c b/src/pntpos.c
index f3a5895..2b006c9 100644
--- a/src/pntpos.c
+++ b/src/pntpos.c
@@ -235,7 +235,8 @@ static int rescode(int iter, const obsd_t *obs, int n, const double *rs,
     for (i=*ns=0;i<n&&i<MAXOBS;i++) {
         vsat[i]=0; azel[i*2]=azel[1+i*2]=resp[i]=0.0;
         
-        if (!(sys=satsys(obs[i].sat,NULL))) continue;
+        if (!(sys=satsys(obs[i].sat,NULL)))
+            continue;
         
         /* reject duplicated observation data */
         if (i<n-1&&i<MAXOBS-1&&obs[i].sat==obs[i+1].sat) {
@@ -246,17 +247,21 @@ static int rescode(int iter, const obsd_t *obs, int n, const double *rs,
         }
         /* geometric distance/azimuth/elevation angle */
         if ((r=geodist(rs+i*6,rr,e))<=0.0||
-            satazel(pos,e,azel+i*2)<opt->elmin) continue;
+            satazel(pos,e,azel+i*2)<opt->elmin)
+            continue;
         
         /* psudorange with code bias correction */
-        if ((P=prange(obs+i,nav,azel+i*2,iter,opt,&vmeas))==0.0) continue;
+        if ((P=prange(obs+i,nav,azel+i*2,iter,opt,&vmeas))==0.0)
+            continue;
         
         /* excluded satellite? */
-        if (satexclude(obs[i].sat,vare[i],svh[i],opt)) continue;
+        if (satexclude(obs[i].sat,vare[i],svh[i],opt))
+            continue;
         
         /* ionospheric corrections */
         if (!ionocorr(obs[i].time,nav,obs[i].sat,pos,azel+i*2,
-                      iter>0?opt->ionoopt:IONOOPT_BRDC,&dion,&vion)) continue;
+                      iter>0?opt->ionoopt:IONOOPT_BRDC,&dion,&vion))
+            continue;
         
         /* GPS-L1 -> L1/B1 */
         if ((lam_L1=nav->lam[obs[i].sat-1][0])>0.0) {
@@ -264,12 +269,21 @@ static int rescode(int iter, const obsd_t *obs, int n, const double *rs,
         }
         /* tropospheric corrections */
         if (!tropcorr(obs[i].time,nav,pos,azel+i*2,
-                      iter>0?opt->tropopt:TROPOPT_SAAS,&dtrp,&vtrp)) {
+                      iter>0?opt->tropopt:TROPOPT_SAAS,&dtrp,&vtrp))
             continue;
-        }
+
         /* pseudorange residual */
         v[nv]=P-(r+dtr-CLIGHT*dts[i*2]+dion+dtrp);
         
+        //printf("sat=%2d:\n",obs[i].sat);
+        //printf("\tx=%5.1f %5.1f %5.1f\n",rr[0],rr[1],rr[2]);
+        //printf("\tnorm=%5.1f\n",r);
+        //printf("\tdtr=%5.8f\n",dtr);
+        //printf("\tP-(dtr-CLIGHT*dts[i*2]+dion+dtrp)=%5.1f\n",P-(dtr-CLIGHT*dts[i*2]+dion+dtrp));
+        //printf("\tP-(-CLIGHT*dts[i*2]+dion+dtrp)=%5.1f\n",P-(-CLIGHT*dts[i*2]+dion+dtrp));
+        //printf("\tsat=%5.1f %5.1f %5.1f\n",*(rs+1*6),*(rs+1*6+1),*(rs+1*6+2));
+        //printf("\terror=%5.1f\n",v[nv]);
+
         /* design matrix */
         for (j=0;j<NX;j++) H[j+nv*NX]=j<3?-e[j]:(j==3?1.0:0.0);
         
@@ -286,6 +300,7 @@ static int rescode(int iter, const obsd_t *obs, int n, const double *rs,
         
         trace(4,"sat=%2d azel=%5.1f %4.1f res=%7.3f sig=%5.3f\n",obs[i].sat,
               azel[i*2]*R2D,azel[1+i*2]*R2D,resp[i],sqrt(var[nv-1]));
+
     }
     /* constraint to avoid rank-deficient */
     for (i=0;i<4;i++) {
-- 
GitLab