- added check for the obs in sat selection.

Author: hirokawa

cssrlib/data/3034265G.21O

@@ -43,11 +43,11 @@ class uGNSS(IntEnum):
     IRN = 7
     GNSSMAX = 8
     GPSMAX = 32
-#    GALMAX = 36
-    GALMAX = 0
+    GALMAX = 36
+    QZSMAX = 10    
+#    GALMAX = 0
 #    BDSMAX = 63
-#    QZSMAX = 10
-    QZSMAX = 0
+#    QZSMAX = 0
 #    GLOMAX = 24
 #    SBSMAX = 24
 #    IRNMAX = 10
@@ -56,8 +56,7 @@ class uGNSS(IntEnum):
     SBSMAX = 0
     IRNMAX = 0
     NONE = -1
-    MAXSAT = GPSMAX+GLOMAX+GALMAX+BDSMAX+QZSMAX
-    # MAXSAT=GPSMAX+GLOMAX+GALMAX+BDSMAX+QZSMAX+SBSMAX+IRNMAX
+    MAXSAT=GPSMAX+GLOMAX+GALMAX+BDSMAX+QZSMAX+SBSMAX+IRNMAX
 
 
 class uSIG(IntEnum):
@@ -177,8 +176,8 @@ def __init__(self):
         self.freq = [1.57542e9, 1.22760e9, 1.17645e9, 1.20714e9]
         self.rb = [0, 0, 0]  # base station position in ECEF [m]
         self.smode = 0  # position mode 0:NONE,1:std,2:DGPS,4:fix,5:float
-        # self.gnss_t = [uGNSS.GPS, uGNSS.GAL, uGNSS.QZS]
-        self.gnss_t = [uGNSS.GPS]
+        self.gnss_t = [uGNSS.GPS, uGNSS.GAL, uGNSS.QZS]
+        # self.gnss_t = [uGNSS.GPS]
         self.loglevel = 2
         self.cnr_min = 35
         self.maxout = 5  # maximum outage [epoch]

cssrlib/data/SEPT2650.21P

@@ -11,7 +11,7 @@
 from cssrlib.ephemeris import satposs
 from cssrlib.cssrlib import sSigGPS, sSigGAL, sSigQZS, sCType
 from cssrlib.ppp import tidedisp, shapiro, windupcorr
-from cssrlib.rtk import IB, ddres, resamb_lambda, valpos, holdamb, initx
+from cssrlib.rtk import IB, ddres, resamb_lambda, valpos, holdamb, initx, kfupdate
 
 MAXITR = 10
 ELMIN = 10
@@ -350,7 +350,7 @@ def ppprtkpos(nav, obs, cs):
     v, H, R = ddres(nav, xp, y, e, sat, el)
     Pp = nav.P.copy()
     # Kalman filter measurement update
-    xp, Pp = gn.kfupdate(xp, Pp, H, v, R)
+    xp, Pp, _ = kfupdate(xp, Pp, H, v, R)
 
     # non-differencial residual for rover after measurement update
     yu, eu, elu = zdres(nav, obs, rs, vs, dts, svh, xp[0:3], cs)

cssrlib/data/SEPT265G.21O

@@ -449,11 +449,35 @@ def udstate(nav, obs, obsb, iu, ir):
 
 def selsat(nav, obs, obsb, elb):
     """ select common satellite between rover and base station """
-    idx0 = np.where(elb >= nav.elmin)
-    idx = np.intersect1d(obs.sat, obsb.sat[idx0], return_indices=True)
+    # exclude satellite with missing observation and cycle slip for rover
+    idx_u = []
+    for k, sat in enumerate(obs.sat):
+        sys, _ = gn.sat2prn(sat)
+        j0 = nav.obs_idx[0][sys]
+        j1 = nav.obs_idx[1][sys]
+        if obs.P[k, j0] == 0.0 or obs.P[k, j1] == 0.0 or \
+            obs.L[k, j0] == 0.0 or obs.L[k, j1] == 0.0 or \
+            obs.lli[k, j0] > 0 or obs.lli[k, j1] > 0:
+            continue
+        idx_u.append(k)
+
+    # exclude satellite with missing observation and cycle slip for base
+    idx_r = []
+    for k, sat in enumerate(obsb.sat):
+        sys, _ = gn.sat2prn(sat)
+        j0 = nav.obs_idx[0][sys]
+        j1 = nav.obs_idx[1][sys]
+        if obsb.P[k, j0] == 0.0 or obsb.P[k, j1] == 0.0 or \
+            obsb.L[k, j0] == 0.0 or obsb.L[k, j1] == 0.0 or \
+            obsb.lli[k, j0] > 0 or obsb.lli[k, j1] > 0 or \
+            elb[k] < nav.elmin:
+            continue
+        idx_r.append(k)
+    
+    idx = np.intersect1d(obs.sat[idx_u], obsb.sat[idx_r], return_indices=True)
     k = len(idx[0])
-    iu = idx[1]
-    ir = idx0[0][idx[2]]
+    iu = np.array(idx_u)[idx[1]]
+    ir = np.array(idx_r)[idx[2]]
     return k, iu, ir
 
 
@@ -498,11 +522,14 @@ def relpos(nav, obs, obsb):
     rsb, _, dtsb, svhb = satposs(obsb, nav)
 
     # non-differencial residual for base
-    yr, er, el = zdres(nav, obsb, rsb, dtsb, svhb, nav.rb, 0)
-    ns, iu, ir = selsat(nav, obs, obsb, el)
+    yr, er, elr = zdres(nav, obsb, rsb, dtsb, svhb, nav.rb, 0)
+    ns, iu, ir = selsat(nav, obs, obsb, elr)
     y = np.zeros((ns*2, nf*2))
     e = np.zeros((ns*2, 3))
 
+    if ns<4:
+        return -1
+    
     y[ns:, :] = yr[ir, :]
     e[ns:, :] = er[ir, :]
 

cssrlib/gnss.py

@@ -13,7 +13,7 @@
 dec = rnxdec()
 nav = stdinit()
 nav = dec.decode_nav(navfile, nav)
-nep = 120
+nep = 360
 t = np.zeros(nep)
 enu = np.zeros((nep, 3))
 sol = np.zeros((nep, nav.nx))

cssrlib/ppprtk.py

@@ -43,6 +43,8 @@
     obs, obsb = rn.sync_obs(dec, decb)
     if ne == 0:
         t0 = nav.t = obs.t
+    if ne>=18:
+        print(ne)
     relpos(nav, obs, obsb)
     t[ne] = gn.timediff(nav.t, t0)
     sol = nav.xa[0:3]

cssrlib/rtk.py

@@ -20,10 +20,9 @@
 decb.decode_obsh(basefile)
 dec.decode_obsh(obsfile)
 
-nep = 85
+nep = 360
 # GSI 3034 fujisawa
 nav.rb = [-3959400.631, 3385704.533, 3667523.111]
-nav.excl_sat = [158]
 t = np.zeros(nep)
 enu = np.zeros((nep, 3))
 smode = np.zeros(nep, dtype=int)
@@ -37,8 +36,6 @@
         obs, obsb = rn.sync_obs(dec, decb)
         if ne == 0:
             t0 = nav.t = obs.t
-        if ne>=60:
-            ne
         relpos(nav, obs, obsb)
         t[ne] = gn.timediff(nav.t, t0)
         sol = nav.x[0:3]