# Patent

## Title
Geofencing and route adherence in global positioning system with signals from fewer than three satellites

## Abstract

Particular implementations are particularly useful in providing a system in which the hardware is more easily upgradable and new hardware functionality may be added without adding any new physical hardware. Through placement of an FPGA closely associated with the CPU of a personal computer, the FPGA may be reconfigured to act as new hardware. A system for installing new virtual hardware involves loading firmware into memory associated with the FPGA and reconfiguring the FPGA through a microcontroller. Particular implementations include universal ports associated with the FPGA into which adapter plugs can be placed to quickly adapt to any device that may be added through the virtual hardware use of the FPGA. Other implementations include high density connectors into which a plurality of ports of varying configurations may be plugged for connection of external electronic equipment through the FPGA.

## Background
FIELD OF THE INVENTION
The present invention is directed to a system, method and article of manufacture for GPS location and more particularly to such a system, method and article of manufacture in which a minimum distance from a GPS receiver to a nearby reference location is determined when there are fewer than three GPS satellites in view.
DESCRIPTION OF RELATED ART
One of the major problems in the use of satellite navigation is blockage of the view of GPS satellites. The blockage is often only partial, with only a few satellites visible and the rest being blocked by trees, buildings and other obstacles. To obtain a GPS fix, a GPS receiver without an on-board atomic clock generally needs to have at least 4 satellites in view (this is called a 3D GPS fix.). With an assumed receiver altitude (typically from a previous 3D fix), only 3 satellites are necessary (so called 2D fix.)
With only two satellites in view, it is impossible to determine the receiver's current location. However, a geometrical surface of the possible user locations, or the Solution Surface, can be determined. In general, knowing the Solution Surface is not very useful, since it is not a bounded area.
SUMMARY OF THE INVENTION
A need therefore exists in the art to provide useful information in situations in which fewer than three satellites are in view.
It is therefore an object of the invention to be able to derive such information.
To achieve the above and other objects, the present invention is directed to a technique in which, for a known reference (base) location (e.g. the receiver's location in the past, or some other position of interest), the distance from this location to the Solution Surface can be computed. That distance is equal to the minimum possible distance between the base location and a current user position. In other words, that distance is how far from the base position the receiver is guaranteed to be.
In a geofence application, the geofence anchor location is set to be the base location; comparing the distance from the anchor location to the Solution Surface with the geofence radius, we can generate a Departure Geofence Crossing Event.
The route adherence algorithm is similar to the geofencing algorithm, except that the minimum-distance test is applied over a set of points or line segments defining a route.
Given a reference (base) geographical location that is known to be relatively close to the current receiver's location compared to the distance to the satellites, the distance from this location to the solution surface can be computed, and that distance is the minimum possible distance from the base location to the current receiver's location. The present disclosure contains the formulae to compute that distance.
In addition to covering the basic techniques for computing the distance from a base point of interest to the solution surface, some geofencing and route adherence applications are discussed and explored. Application of those algorithms would provide better availability of geofence monitoring and route adherence testing in limited GPS satellite visibility environments.
In other embodiments, an alternate derivation for the linearized solution space equation is shown, and a similar positioning technique with a single satellite and a high accuracy terminal clock is detailed.
Throughout the present disclosure, the following terminology will be used:
User and Receiver are used interchangeably.
Location and Position are used interchangeably.
ECEF—Earth centered, Earth fixed coordinate system (e.g. WGS-84 or ITRF2000)
Base Location—A reference location the distance to which we are interested in.
Geofence Anchor Position—The geographical position where geofence monitored object is expected to be on average.
User Coordinate Frame—a coordinate system rotated so that the basis vectors point in cardinal directions (e.g. east, north and up) at the user or base location. The origin of the coordinate system does not need to be at the user or base location.
Locally Horizontal Plane—A horizontal plane in the user coordinate frame at the base location. This plane is parallel to the WGS-84 ellipsoid at the user location (tangent if the altitude is zero with respect to the WGS-84 ellipsoid).
Solution Surface—A geometrical surface of possible user positions obtained from the GPS equations for two satellites.
Solution Plane—A linearized approximation to the solution surface.
Horizontal Solution Line—The line formed by the intersection of the solution plane with the locally horizontal plane.
Departure Geofence Crossing Event (DGCE)—The event generated by a geofence software application indicating that a tracked object is moving outside of the circle with a given geofence radius.
In addition, throughout the present disclosure, the following notation will be used:
{right arrow over (x)}S,i=[xS,i yS,i zS,i] The i-th satellite position vector (typically but not necessarily in ECEF coordinates).
{right arrow over (x)}U=[xU yU zU] The user's position vector (also typically in ECEF coordinates).
{right arrow over (r)}i={right arrow over (x)}S,i−{right arrow over (x)}U The range vector from the user position to the i-th satellite.
{right arrow over (l)}i={right arrow over (r)}i/|{right arrow over (r)}i| The line of sight vector from the user position to the i-th satellite position.
{right arrow over (x)}B The base position. This is often a previous user position but could be any location of interest in general.
{right arrow over (r)}i,B The range vector from base position to the i-th satellite.
{right arrow over (l)}B The line-of-sight vector from the base position to the i-th satellite position.
{right arrow over (δ)}U={right arrow over (x)}U−{right arrow over (x)}B The displacement vector from the base position to the current user position.
bU The clock bias of the receiver's clock.
ρi=ri+bU The pseudorange to the i-th satellite.
NS The number of satellites used in computation.
D(SPl,B) The distance from the base position to the Solution Plane.
D(SLn,B) The distance from the base position to the Solution Line.
λS1,2=∠{right arrow over (l)}1{right arrow over (l)}2 The angle between the line-of-sight vectors of SV1 and SV2.
λS1,2min The minimum angle between SV1 and SV2 allowing geofence monitoring.
Fconf The confidence factor used in the computation of λS1,2min.
σρ The root-mean-square error of a pseudorange estimate.
σρ1,2=√{square root over (2)}σρ The root-mean-square error of the pseudorange difference of ρ2−ρ1.
{right arrow over (n)}LHPl,B The unitary vector normal to the locally horizontal plane at the base location.
Δρ=ρ2−ρ1 The pseudorange difference.
ΔrB=r2,B−r1,B The true range difference at the base location.
λ(SPl,LHPl) The dihedral angle between the Solution Plane and the Locally Horizontal Plane.



BRIEF DESCRIPTION OF THE DRAWINGS
Preferred embodiments of the invention will be set forth in detail with reference to the drawings, in which:
FIGS. 1A and 1B are diagrams which show examples of a Solution Surface for a first preferred embodiment;
FIG. 2 is a diagram which shows the derivation of a solution line equation in a case in which the location of the receiver is assumed to be horizontally constrained;
FIG. 3 is a flow chart showing the computation of the minimum distance and its comparison to a given geofence radius;
FIG. 4 is a flow chart showing the computation of the minimum distance to one of a plurality of line segments in a route;
FIG. 5 is a block diagram showing an implementation of any of the preferred embodiments; and
FIG. 6 is a diagram which shows an example of a Solution Surface for another preferred embodiment.



DETAILED DESCRIPTION OF THE PREFERRED EMBODIMENTS
Preferred embodiments of the invention will now be set forth in detail with reference to the drawings, in which like reference numerals refer to like elements throughout.
In a first preferred embodiment, the distance is determined to a user position on a plane in the case of two satellites and a GPS receiver with a standard clock. The basic GPS fix equation can be written as:

{ρi=ri+bU|i=1, 2, . . . NS}  (E10)

where

ri=|{right arrow over (x)}S,i−{right arrow over (x)}U|=√{square root over (({right arrow over (x)}S,i−{right arrow over (x)}U)·({right arrow over (x)}S,i−{right arrow over (x)}U))}  (E20)

In Cartesian coordinates the equations in (E10) will take the following form:




  
    
      
        
          ρ
          i
        
        =
        
          
            
              
                
                  (
                  
                    
                      x
                      
                        S
                        ,
                        i
                      
                    
                    -
                    
                      x
                      U
                    
                  
                  )
                
                2
              
              +
              
                
                  (
                  
                    
                      y
                      
                        S
                        ,
                        i
                      
                    
                    -
                    
                      y
                      U
                    
                  
                  )
                
                2
              
              +
              
                
                  (
                  
                    
                      z
                      
                        S
                        ,
                        i
                      
                    
                    -
                    
                      z
                      U
                    
                  
                  )
                
                2
              
            
          
          +
          
            b
            U
          
        
      
    
    
      
        (
        E30
        )
      
    
  




Here ρi are known measurements, {right arrow over (x)}S,i are computed from GPS navigation and observation data values, and {right arrow over (x)}U and bU are unknowns.
With two satellites in view, equation (E30) is a pair of dependent equations. By equating clock biases bU in (E30) for the two satellites we obtain:




  
    
      
        
          
            ρ
            1
          
          -
          
            
              
                
                  (
                  
                    
                      x
                      
                        S
                        ,
                        1
                      
                    
                    -
                    
                      x
                      U
                    
                  
                  )
                
                2
              
              +
              
                
                  (
                  
                    
                      y
                      
                        S
                        ,
                        1
                      
                    
                    -
                    
                      y
                      U
                    
                  
                  )
                
                2
              
              +
              
                
                  (
                  
                    
                      z
                      
                        S
                        ,
                        1
                      
                    
                    -
                    
                      z
                      U
                    
                  
                  )
                
                2
              
            
          
        
        =
        
          
            ρ
            2
          
          -
          
            
              
                
                  (
                  
                    
                      x
                      
                        S
                        ,
                        2
                      
                    
                    -
                    
                      x
                      U
                    
                  
                  )
                
                2
              
              +
              
                
                  (
                  
                    
                      y
                      
                        S
                        ,
                        2
                      
                    
                    -
                    
                      y
                      U
                    
                  
                  )
                
                2
              
              +
              
                
                  (
                  
                    
                      z
                      
                        S
                        ,
                        1
                      
                      
                        (
                        2
                        )
                      
                    
                    -
                    
                      z
                      U
                    
                  
                  )
                
                2
              
            
          
        
      
    
    
      
        (
        E40
        )
      
    
  




Equation (E40) thus defines the Solution Surface.
The two dimensional profile of this surface (intersection with a plane) is a curve. In a coordinate system where both satellites and the receiver lie in the XY plane we can write the profile curve equation as:




  
    
      
        
          
            ρ
            1
          
          -
          
            
              
                
                  (
                  
                    
                      x
                      
                        S
                        ,
                        1
                      
                    
                    -
                    
                      x
                      U
                    
                  
                  )
                
                2
              
              +
              
                
                  (
                  
                    
                      y
                      
                        S
                        ,
                        1
                      
                    
                    -
                    
                      y
                      U
                    
                  
                  )
                
                2
              
            
          
        
        =
        
          
            ρ
            2
          
          -
          
            
              
                
                  (
                  
                    
                      x
                      
                        S
                        ,
                        2
                      
                    
                    -
                    
                      x
                      U
                    
                  
                  )
                
                2
              
              +
              
                
                  (
                  
                    
                      y
                      
                        S
                        ,
                        2
                      
                    
                    -
                    
                      y
                      U
                    
                  
                  )
                
                2
              
            
          
        
      
    
    
      
        (
        E50
        )
      
    
  




We therefore have an implicit function yu=ƒ(xu) that can be plotted to gain some intuitive understanding of the Solution Surface's behaviour.
FIGS. 1A and 1B show examples of the Solution Surface profile curves for this case. In both of those figures, there are only two satellites visible, namely SV1 and SV2. The actual user's position AUP on the earth's surface ES cannot be determined unequivocally. However, a solution surface SS can be determined as disclosed above. Further, if the actual user's position AUP is known to be on the earth's surface ES, as opposed to being at some other location such as an airplane over the earth's surface, the intersection of the solution surface SS and the earth's surface ES is a curve on which the actual user's position AUP must be located.
To derive the solution plane equation, we linearize equations (E10) around the base location {right arrow over (x)}B.




  
    r
    i
  
  =
  
    
      
      
        
          
            x
            →
          
          
            S
            ,
            i
          
        
        -
        
          
            x
            →
          
          U
        
      
      
    
    =
    
      
        
        
          
            
              x
              →
            
            
              S
              ,
              i
            
          
          -
          
            
              x
              →
            
            B
          
          +
          
            
              δ
              →
            
            U
          
        
        
      
      =
      
        
          
          
            
              
                r
                →
              
              
                i
                ,
                B
              
            
            +
            
              δ
              U
            
          
          
        
        =
        
          
            
              
                r
                
                  i
                  ,
                  B
                
                2
              
              +
              
                2
                ⁢
                
                  r
                  
                    i
                    ,
                    B
                  
                
                ⁢
                
                  
                    
                      δ
                      →
                    
                    U
                  
                  ·
                  
                    
                      l
                      →
                    
                    
                      i
                      ,
                      B
                    
                  
                
              
              +
              
                δ
                U
                2
              
            
          
          =
          
            
              r
              
                i
                ,
                B
              
              
                (
                i
                )
              
            
            ⁢
            
              
                1
                +
                
                  2
                  ⁢
                  
                    
                      
                        
                          δ
                          →
                        
                        U
                      
                      
                        r
                        
                          i
                          ,
                          B
                        
                      
                    
                    ·
                    
                      
                        l
                        →
                      
                      
                        i
                        ,
                        B
                      
                    
                  
                
                +
                
                  
                    (
                    
                      
                        δ
                        U
                      
                      
                        r
                        
                          i
                          ,
                          B
                        
                      
                    
                    )
                  
                  2
                
              
            
          
        
      
    
  




For linearization, the distance to a satellite is assumed to be much greater than the displacement from the base location:




  
    
      
      
        δ
        U
      
      
    
    
      
      
        r
        
          i
          ,
          B
        
      
      
    
  
  ⁢
  
    <<
    1
  




We can then use the Taylor series approximation for the square root function:

√{square root over (1+α)}≈1+1/2α; |α|<<1

so

ri≈ri,B+{right arrow over (δ)}U·{right arrow over (l)}i,B

And

ρi≈ri,B+{right arrow over (δ)}U·{right arrow over (l)}i,B+bU  (E60)

In the two-satellite case this is a pair of equations. Eliminating the clock bias by equating it in the two equations for SV1 and SV2 in (E60) we obtain:

ρ1−r1=ρ2−r2

ρ1−(r1,B+{right arrow over (δ)}U·{right arrow over (l)}1,B)=ρ2−(r2,B+{right arrow over (δ)}U·{right arrow over (l)}2,B)

{right arrow over (δ)}U·({right arrow over (l)}2,B−{right arrow over (l)}1,B)={right arrow over (δ)}U·Δ{right arrow over (l)}B=(ρ2−ρ1)−(r2,B−r1,B)=Δρ−ΔrB

Dividing both parts by |Δ{right arrow over (l)}B| we have:




  
    
      
        
          
            
              
                δ
                →
              
              U
            
            ·
            
              
                Δ
                ⁢
                
                  
                
                ⁢
                
                  
                    l
                    →
                  
                  B
                
              
              
                
                
                  Δ
                  ⁢
                  
                    
                  
                  ⁢
                  
                    l
                    B
                  
                
                
              
            
          
          =
          
            
              
                Δρ
                -
                
                  Δ
                  ⁢
                  
                    
                  
                  ⁢
                  
                    r
                    B
                  
                
              
              
                
                
                  Δ
                  ⁢
                  
                    
                  
                  ⁢
                  
                    l
                    B
                  
                
                
              
            
            =
            
              -
              
                D
                
                  (
                  
                    SPI
                    ,
                    B
                  
                  )
                
              
            
          
        
        ⁢
        
          

        
        ⁢
        Or
      
    
    
      
        (
        E70
        )
      
    
  
  
    
      
        
          
            (
            
              
                
                  x
                  →
                
                U
              
              -
              
                
                  x
                  →
                
                B
              
            
            )
          
          ·
          
            
              Δ
              ⁢
              
                
              
              ⁢
              
                
                  l
                  →
                
                B
              
            
            
              
              
                Δ
                ⁢
                
                  
                
                ⁢
                
                  
                    l
                    →
                  
                  B
                
              
              
            
          
        
        =
        
          -
          
            D
            
              (
              
                SPl
                ,
                B
              
              )
            
          
        
      
    
    
      
        (
        E71
        )
      
    
  




The obtained equation (E71) is the Hessian normal form of the plane equation, where |D(SPl,B)| is the distance from the plane of possible user positions {right arrow over (x)}u to the point {right arrow over (x)}B·|D(SPl,B)| is the minimum possible distance from the base position {right arrow over (x)}B to the unknown current user position {right arrow over (x)}u.
We can rewrite D(SPl,B):




  
    
      
        
          D
          
            (
            
              SPl
              ,
              B
            
            )
          
        
        =
        
          
            
              
                Δ
                ⁢
                
                  
                
                ⁢
                
                  r
                  B
                
              
              -
              Δρ
            
            
              
              
                Δ
                ⁢
                
                  
                
                ⁢
                
                  
                    l
                    →
                  
                  B
                
              
              
            
          
          =
          
            
              
                
                  Δ
                  ⁢
                  
                    
                  
                  ⁢
                  
                    r
                    B
                  
                
                -
                Δρ
              
              
                2
                ·
                
                  (
                  
                    1
                    -
                    
                      cos
                      ⁢
                      
                        
                      
                      ⁢
                      
                        ∠
                        ⁡
                        
                          (
                          
                            
                              
                                l
                                →
                              
                              
                                1
                                ,
                                B
                              
                            
                            ,
                            
                              
                                l
                                →
                              
                              
                                2
                                ,
                                B
                              
                            
                          
                          )
                        
                      
                    
                  
                  )
                
              
            
            =
            
              
                
                  Δ
                  ⁢
                  
                    
                  
                  ⁢
                  
                    r
                    B
                  
                
                -
                Δρ
              
              
                2
                ·
                
                  (
                  
                    1
                    -
                    
                      cos
                      ⁢
                      
                        
                      
                      ⁢
                      
                        λ
                        
                          
                            S
                            
                              1
                              ,
                              2
                            
                          
                          ,
                          B
                        
                      
                    
                  
                  )
                
              
            
          
        
      
    
    
      
        (
        E80
        )
      
    
  




A solution line equation will now be disclosed. In many cases it may be useful to assume that the user positions are constrained within the locally horizontal plane defined at the same altitude as the base position. In cases where the Solution Plane is significantly inclined with respect to the locally horizontal plane, and where motion is essentially horizontal, the minimum distance from the base position to the Solution Plane will greatly underestimate the actual distance of the user to the Solution Plane.
With the assumption that the receiver is constrained to move within the Locally Horizontal Plane, as shown in FIG. 2, we have the receiver's solution space narrowed to the line (Solution Line SL) formed by the intersection of the Solution Plane SP and the Locally Horizontal Plane LHP. The distance to the Solution Line SL can be computed from the diagram of FIG. 2 by the following equation:




  
    
      
        
          D
          
            (
            
              SLn
              ,
              B
            
            )
          
        
        =
        
          
            
              D
              
                (
                
                  SPl
                  ,
                  B
                
                )
              
            
            
              sin
              ⁢
              
                
              
              ⁢
              ∠
              ⁢
              
                
              
              ⁢
              Δ
              ⁢
              
                
              
              ⁢
              
                
                  l
                  →
                
                B
              
              ⁢
              
                
                  n
                  →
                
                
                  LHPl
                  ,
                  B
                
              
            
          
          =
          
            
              
                Δ
                ⁢
                
                  
                
                ⁢
                
                  r
                  B
                
              
              -
              
                Δ
                ⁢
                
                  
                
                ⁢
                ρ
              
            
            
              
                2
                ·
                
                  (
                  
                    1
                    -
                    
                      cos
                      ⁢
                      
                        
                      
                      ⁢
                      
                        λ
                        
                          
                            SV
                            
                              1
                              ,
                              2
                            
                          
                          ,
                          B
                        
                      
                    
                  
                  )
                
                ·
                sin
              
              ⁢
              
                
              
              ⁢
              
                λ
                
                  (
                  
                    SPl
                    ,
                    LHPl
                  
                  )
                
              
            
          
        
      
    
    
      
        (
        E90
        )
      
    
  




Note: the receiver position's constraining plane does not necessarily have to be the locally horizontal plane; however, this is the most likely scenario in realistic applications.
The sensitivity of D(SPl,B) and D(SLn,B) to pseudorange errors will now be considered. In satellite navigation it is common to define a constellation geometry quality factor known as dilution of precision. This unitless quantity is multiplied by the expected root-mean-square pseudorange error (assumed to be the same for all satellites) to estimate the positioning error.
The equivalents of the position dilution of precision (PDOP) used in a regular GPS fix in our case would be (for the distances to the Solution Plane and Solution Line respectively):




  
    
      
        
          Γ
          Pl
        
        =
        
          
            2
          
          
            2
            ·
            
              (
              
                1
                -
                
                  cos
                  ⁢
                  
                    
                  
                  ⁢
                  
                    λ
                    
                      
                        SV
                        
                          1
                          ,
                          2
                        
                      
                      ,
                      B
                    
                  
                
              
              )
            
          
        
      
    
    
      
        (
        E91
        )
      
    
  
  
    
      
        
          Γ
          Ln
        
        =
        
          
            2
          
          
            
              2
              ·
              
                (
                
                  1
                  -
                  
                    cos
                    ⁢
                    
                      
                    
                    ⁢
                    
                      λ
                      
                        
                          SV
                          
                            1
                            ,
                            2
                          
                        
                        ,
                        B
                      
                    
                  
                
                )
              
              ·
              sin
            
            ⁢
            
              
            
            ⁢
            
              λ
              
                (
                
                  SPl
                  ,
                  LHPl
                
                )
              
            
          
        
      
    
    
      
        (
        E92
        )
      
    
  




The pseudorange estimate mean error would translate in the mean D(SPl,B) error as:

σD(SPl,B)=ΓPl·σρ  (E93)

Similarly, the pseudorange estimate mean error would translate in the mean DSLn,B error as:

σD(SLn,B)=ΓLn·σρ  (E93)

The √{square root over (2)} value in the numerators of E91 and E92 appears from an assumption that the pseudorange errors from each satellite are uncorrelated and have the same root-mean-square value σρ, in which case

σρ1,2=√{square root over (2)}·σρ  (E95)

A table of sample values for ΓPl(λSV1,2) is provided below, showing its dependence on the angle λSV1,2 between the two satellite line-of-sight vectors.







TABLE 1







Values of Γ(λSV1,2)









λSV1,2
















1°
5°
7°
10°
20°
30°
45°
90°



















Γ
4642.71
185.82
94.86
46.54
11.73
5.28
2.41
0.707









The table shows that the accuracy decreases dramatically for low values of λSV1,2.
The departure geofence algorithm will now be derived. The purpose of the proposed algorithm is to generate an alarm when a tracked object departs from a known geofence anchor location (base location) for more than an allowed distance (geofence radius). The algorithm assumes only 2 GPS satellites are in view.
The algorithm is based on the computation of the D(SPl,B) value and comparing its absolute value against a given geofence radius Rgeof. The algorithm comprises the following steps, shown in the flow chart of FIG. 3:
Step 302: Decide whether we can use the current satellite geometry. Compute λSV1,2,B. If λSV1,2,B<λSV1,2min then D(SPl,B) shall be assumed incomputable in step 304, and geofence monitoring will be declared unavailable.
λSV1,2min is a function of Rgeof:




  
    
      
        
          
            
              σ
              ρ
            
            
              
                2
              
              ·
              
                (
                
                  1
                  -
                  
                    cos
                    ⁢
                    
                      
                    
                    ⁢
                    
                      λ
                      
                        SV
                        
                          1
                          ,
                          2
                        
                      
                      min
                    
                  
                
                )
              
            
          
          =
          
            
              R
              geof
            
            ⁢
            
              F
              conf
            
          
        
        ⁢
        
          

        
        ⁢
        Hence
        ⁢
        
          

        
        ⁢
        
          
            λ
            
              SV
              
                1
                ,
                2
              
            
            min
          
          =
          
            arccos
            (
            
              1
              -
              
                
                  σ
                  ρ
                
                
                  
                    2
                  
                  ⁢
                  
                    R
                    geof
                  
                  ⁢
                  
                    F
                    conf
                  
                
              
            
            )
          
        
      
    
    
      
        (
        E100
        )
      
    
  




Step 306: Compute D(SPl,B). If it is determined in step 308 that

|D(SPl,B)|>Rgeof  (E110)

generate the departure geofence crossing event in step 310. Otherwise, in step 312, the geofence crossing event is not generated.
The algorithm can be modified for use of D(SLn,B) instead of D(SPl,B).
The performance of the algorithm will now be discussed. The probability of a missed geofence crossing event is quite large due to a ‘stealth’ factor. When the receiver is moving within the Solution Plane, the computed D(SPl,B) value does not change and the receiver's movement remains undetected. The ‘stealth’ factor contribution to the missed detection rate is proportional to




  1
  
    R
    geof
    2
  




for the receiver moving on the horizontal plane. The other factors contributing to missed detection are the pseudorange and satellite position errors.

The probability of false geofence crossing events (false alarms) is fully determined by the errors in D(SPl,B) caused by the errors in pseudorange estimates and satellite positions.
In a second preferred embodiment, a route adherence algorithm is used. The previous algorithm for geofencing, in which the minimum possible distance from a reference (base) position is computed, can be extended to route adherence testing by computing the minimum possible distance for a number of points along the route.
One approach is to approximate the route to be adhered to using a number of line segments in three-dimensional space, where the lengths of the line segments are short compared to the distance to the satellites (i.e. less than a few hundred kilometers). In that case, the linearized solution planes at each end of the line segment are approximately identical. For each line segment, there are four possible cases:
1. The line segment passes through the solution plane;
2. The line segment does not pass through the solution plane, and the first point on the segment is closest to the solution plane; or
3. The line segment does not pass through the solution plane, and the first point on the segment is closest to the solution plane; or
4. The line segment is exactly parallel to the solution plane
In practice, the last case is extremely unlikely to occur, and for the purposes of determining the minimum distance over the segment, is equivalent to either case 2 or 3.
The route adherence algorithm computes the closest distance of approach of any solution plane to the entire route. If the minimum possible distance is greater than a selected threshold distance, then it is known that the user has not adhered to the route.
The minimum distance algorithm for a single line segment is as follows, as shown in the flow chart of FIG. 4:
Step 402: For the first point in each line segment, compute d0=D(SPl,B) and the normal vector




  
    
      
        v
        →
      
      0
    
    =
    
      
        Δ
        ⁢
        
          
        
        ⁢
        
          
            l
            →
          
          B
        
      
      
        sign
        ⁢
        
          
        
        ⁢
        
          D
          
            (
            
              SPl
              ,
              B
            
            )
          
        
      
    
  
  ,




with the first point of the segment used as the base

Step 404: For the second point in each line segment, compute d1=D(SPl,B) and the normal vector




  
    
      
        v
        →
      
      1
    
    =
    
      
        Δ
        ⁢
        
          
        
        ⁢
        
          
            l
            →
          
          B
        
      
      
        sign
        ⁢
        
          
        
        ⁢
        
          D
          
            (
            
              SPl
              ,
              B
            
            )
          
        
      
    
  
  ,




with the second point of the segment used as the base

Step 406: Determine whether the line segment passes through the solution plane by computing the sign of {right arrow over (ν)}0·{right arrow over (ν)}1. If it is positive, the line segment doesn't pass through the solution plane. If it is negative, it is determined in step 408 that the line segment passes through the solution plane, and the minimum distance over the line segment is 0.
Step 410: If the line segment doesn't pass through the solution plane, then the minimum distance of the line segment to the solution plane is min(d0, d1).
This basic algorithm could be applied over the full set of line segments describing the route in a number of ways. For example, the algorithm could be applied over each line segment and the global minimum could be found. Alternatively, the set of line segments to be searched could be restricted to a subset using side information such as that last known user position and the maximum possible speed.
No hardware modifications are necessary to realize the present invention in an existing GPS receiver. As for software modifications, the GPS fix computation software module should be augmented with the computation of the defined above λSV1,2 and D(SPl,B) when only two satellites are in view. As shown in FIG. 5, a GPS receiver 500 includes an antenna 502, a microprocessor 504 and a storage medium 506, in communication with or integrated into the microprocessor 504, for storing the software. The software can be supplied to the receiver 500 by way of any suitable medium 508. Geofence API functions utilizing the above values should be provided to the geofence application.
As for which parameter to use, D(SPl,B) or D(SLn,B), the user should use the horizontal case if the tracked object is expected to remain in the horizontal plane. The use of D(SPl,B) is more generic.
The situation where only two satellites are visible quite often will correlate with a situation in which only a small aperture of the sky is visible; therefore, the angle between the two visible satellites will often be small. The disadvantage of the small angle between the two satellites is partially compensated by smaller ionospheric and tropospheric errors due to the fact that the ionospheric and tropospheric delays for both satellites will be close.
In the case of an unknown base location, instead of using the base location, we can use a previously computed Solution Plane (base Solution Plane) and measure the minimum distance between the base and current planes.
The best accuracy (precision) would be achieved if the base Solution Plane is originated by the same SVs as the current Solution Plane and the time difference between the computations is not more than a couple of minutes.
Pseudorange dynamic filtering for static receivers can be beneficial in removing thermal noise errors, especially in the case of small angle λSV1,2 where ionospheric error contribution is small. It probably will also filter out free motion receiver's movement around the anchor position.
Another approach to finding the solution plane is to approximate pseudoranges for the current user position with the pseudoranges for the base position using the formulae for the differential of a multivariable function:




  
    
      ρ
      i
    
    ≈
    
      
        ρ
        
          i
          ,
          B
        
      
      +
      
        
          
            (
            
              
                ∂
                
                  ∂
                  x
                
              
              ⁢
              
                r
                i
              
            
            )
          
          B
        
        ⁢
        
          (
          
            
              x
              U
            
            -
            
              x
              B
            
          
          )
        
      
      +
      
        
          
            (
            
              
                ∂
                
                  ∂
                  y
                
              
              ⁢
              
                r
                i
              
            
            )
          
          B
        
        ⁢
        
          (
          
            
              y
              U
            
            -
            
              y
              B
            
          
          )
        
      
      +
      
        
          
            (
            
              
                ∂
                
                  ∂
                  z
                
              
              ⁢
              
                r
                i
              
            
            )
          
          B
        
        ·
        
          (
          
            
              z
              U
            
            -
            
              z
              B
            
          
          )
        
      
      +
      
        
          
            (
            
              
                ∂
                
                  ∂
                  b
                
              
              ⁢
              
                ρ
                i
              
            
            )
          
          B
        
        ·
        
          (
          
            
              b
              U
            
            -
            
              b
              
                U
                ,
                B
              
            
          
          )
        
      
    
  
  =
  
    
      
        r
        
          i
          ,
          B
        
      
      +
      
        b
        
          U
          ,
          B
        
      
      +
      
        
          
            
              x
              
                S
                ,
                i
              
            
            -
            
              x
              B
            
          
          
            r
            
              i
              ,
              B
            
          
        
        ·
        
          (
          
            
              x
              U
            
            -
            
              x
              B
            
          
          )
        
      
      +
      
        
          
            
              y
              
                S
                ,
                i
              
            
            -
            
              y
              B
            
          
          
            r
            
              i
              ,
              B
            
          
        
        ·
        
          (
          
            
              y
              U
            
            -
            
              y
              B
            
          
          )
        
      
      +
      
        
          
            
              z
              
                S
                ,
                i
              
            
            -
            
              z
              B
            
          
          
            r
            
              i
              ,
              B
            
          
        
        ·
        
          (
          
            
              z
              U
            
            -
            
              z
              B
            
          
          )
        
      
      +
      
        b
        U
      
      -
      
        b
        
          U
          ,
          B
        
      
    
    =
    
      
        r
        
          i
          ,
          B
        
      
      +
      
        
          
            δ
            →
          
          U
        
        ·
        
          
            i
            →
          
          
            i
            ,
            B
          
        
      
      +
      
        b
        U
      
    
  




The obtained equation is the same as (E60).
The case of a GPS receiver with a normal accuracy on-board clock and two GPS satellites in view should be similar to the case of a GPS receiver with a high accuracy ideal (e.g. atomic) clock and a single GPS satellite in view. In the latter case, one satellite in view is sufficient to determine a distance from the base location to the Solution Surface. Since the receiver's clock ideally does not have a clock bias, we can estimate the true range to the satellite, rather than the pseudorange. The difference of the range to the current position and the range to the base location is the distance from the base location to the Solution Surface, which happens to be a sphere:

D(SS,B)=|rcur−rB|

If the base location is the receiver's previous location, D(SS,B) has the meaning of the minimum displacement of the receiver from the previous location.
FIG. 6 shows the distance from the base location BL to the Solution Surface SS (which is a sphere) in the case of one satellite SV and a receiver with an ideal clock. The Solution Surface SS is the set of all possible user locations PUL. Any of the variations disclosed above for the preferred embodiment using two satellites can be implemented as well with the preferred embodiment using only one satellite.
While various preferred embodiments of the invention have been disclosed in detail above, those skilled in the art who have reviewed the present disclosure will readily appreciate that other embodiments can be realized within the scope of the invention. For example, a GPS receiver according to the present invention can either have a simple software upgrade from existing receivers or be purpose built. Also, the present invention is not limited to a single example of a global positioning system, but can instead be used with other global positioning systems, such as GLONASS or the Galileo positioning system. Therefore, the present invention should be construed as limited only by the appended claims.

## Claims

1. A method for determining a minimum distance from a global positioning system receiver to a reference location, the method comprising:
(a) receiving at least one global positioning system signal from no more than two global positioning system satellites in the global positioning system receiver;
(b) determining, from the at least one global positioning system signal received in step (a), a set of possible locations at which the global positioning system can be located; and
(c) determining a minimum distance from the set of possible locations to the reference location as the minimum distance from the global positioning system receiver to the reference location
wherein, in step (a), said at least one global positioning system signal is received from two said global positioning system satellites and wherein step (b) comprises determining a solution plane approximating a solution surface comprising the set of possible locations near the reference location.

2. The method of claim 1, wherein step (c) comprises determining a minimum distance from the solution plane to the reference location.

3. The method of claim 1, wherein the global positioning system receiver is assumed to have a location constraint, and wherein the set of possible locations is determined from the solution plane and the location constraint.

4. The method of claim 3, wherein the location constraint is that the global positioning system receiver is assumed to be on a constraint surface, and wherein the set of possible locations is determined from an intersection of the constraint surface and the solution plane.

5. The method of claim 4, wherein the constraint surface is a horizontal plane.

6. A method for determining a minimum distance from a global positioning system receiver to a reference location, the method comprising:
(a) receiving at least one global positioning system signal from no more than two global positioning system satellites in the global positioning system receiver;
(b) determining, from the at least one global positioning system signal received in step (a), a set of possible locations at which the global positioning system can be located; and
(c) determining a minimum distance from the set of possible locations to the reference location as the minimum distance from the global positioning system receiver to the reference location, wherein the reference location is determined from either a previously determined location of the global positioning system receiver or a previously determined solution surface.

7. A method for determining a minimum distance from a global positioning system receiver to a reference location, the method comprising:
(a) receiving at least one global positioning system signal from no more than two global positioning system satellites in the global positioning system receiver;
(b) determining, from the at least one global positioning system signal received in step (a), a set of possible locations at which the global positioning system can be located; and
(c) determining a minimum distance from the set of possible locations to the reference location as the minimum distance from the global positioning system receiver to the reference location, wherein the reference location is a geofence anchor location.

8. The method of claim 7, wherein the geofence anchor location has a corresponding geofence radius, and further comprising comparing the minimum distance to the geofence radius to determine whether a geofence crossing even has taken place.

9. A method for determining a minimum distance from a global positioning system receiver to a reference location, the method comprising:
(a) receiving at least one global positioning system signal from no more than two global positioning system satellites in the global positioning system receiver;
(b) determining, from the at least one global positioning system signal received in step (a), a set of possible locations at which the global positioning system can be located; and
(c) determining a minimum distance from the set of possible locations to the reference location as the minimum distance from the global positioning system receiver to the reference location, wherein the reference location is on a route, and further comprising using the minimum distance to determine whether the global positioning system receiver is following the route.

10. The method of claim 9, wherein the route is modeled as a plurality of line segments, and wherein the minimum distance is calculated to each of the line segments.

11. A global positioning system receiver for determining a minimum distance to a reference location, the global positioning system receiver comprising:
an antenna for receiving at least one global positioning system signal from no more than two global positioning system satellites; and
a processor, in communication with the antenna, for determining, from the at least one global positioning system signal received by the antenna, a set of possible locations at which the global positioning system can be located; and for determining a minimum distance from the set of possible locations to the reference location as the minimum distance from the global positioning system receiver to the reference location, wherein said at least one global positioning system signal is received from two said global positioning system satellites and wherein the processor determines a solution plane approximating a solution surface comprising the set of possible locations near the reference location.

12. The global positioning system receiver of claim 11, wherein the processor determines a minimum distance from the solution plane to the reference location.

13. The global positioning system receiver of claim 11, wherein the global positioning system receiver is assumed to have a location constraint, and wherein the set of possible locations is determined from the solution plane and the location constraint.

14. The global positioning system receiver of claim 13, wherein the location constraint is that the global positioning system receiver is assumed to be on a constraint surface, and wherein the set of possible locations is determined from an intersection of the constraint surface and the solution plane.

15. The global positioning system receiver of claim 14, wherein the constraint surface is a horizontal plane.

16. A global positioning system receiver for determining a minimum distance to a reference location, the global positioning system receiver comprising:
an antenna for receiving at least one global positioning system signal from no more than two global positioning system satellites; and
a processor, in communication with the antenna, for determining, from the at least one global positioning system signal received by the antenna, a set of possible locations at which the global positioning system can be located; and for determining a minimum distance from the set of possible locations to the reference location as the minimum distance from the global positioning system receiver to the reference location, wherein the reference location is either a previously determined location of the global positioning system receiver or a previously determined solution surface.

17. A global positioning system receiver for determining a minimum distance to a reference location, the global positioning system receiver comprising:
an antenna for receiving at least one global positioning system signal from no more than two global positioning system satellites; and
a processor, in communication with the antenna, for determining, from the at least one global positioning system signal received by the antenna, a set of possible locations at which the global positioning system can be located; and for determining a minimum distance from the set of possible locations to the reference location as the minimum distance from the global positioning system receiver to the reference location, wherein the reference location is a geofence anchor location.

18. The global positioning system receiver of claim 17, wherein the geofence anchor location has a corresponding geofence radius, and further comprising comparing the minimum distance to the geofence radius to determine whether a geofence crossing even has taken place.

19. A global positioning system receiver for determining a minimum distance to a reference location, the global positioning system receiver comprising:
an antenna for receiving at least one global positioning system signal from no more than two global positioning system satellites; and
a processor, in communication with the antenna, for determining, from the at least one global positioning system signal received by the antenna, a set of possible locations at which the global positioning system can be located; and for determining a minimum distance from the set of possible locations to the reference location as the minimum distance from the global positioning system receiver to the reference location, wherein the reference location is on a route, and further comprising using the minimum distance to determine whether the global positioning system receiver is following the route.

20. The global positioning system receiver of claim 19, wherein the route is modeled as a plurality of line segments, and wherein the minimum distance is calculated to each of the line segments.