Resilient Distributed Positioning Networks

Abstract

Co-channel beacon transmissions are provided with at least one of spectral redundancy and temporal redundancy. A receiver produces a snapshot of a superposition of received co-channel beacon transmissions. Subcarrier demodulation, code nulling, or a Class-C linear minimum-mean-square error (MMSE) operation separates multiples ones of the co-channel beacon transmissions or eliminates inter-symbol interference and inter-subcarrier interference in the snapshot. Receiver operations can be performed at a network user, a network node, or a network operations center.

Claims

1. A method performed at a receiver, wherein the receiver is configured to receive a superposition of co-channel beacon transmissions transmitted from a plurality of beacon transmitters, each of the co-channel beacon transmissions having at least one of spectral redundancy and temporal redundancy; the method comprising: generating a snapshot of the superposition; and using subcarrier demodulation, code nulling, or a Class-C linear minimum-mean-square error (MMSE) operation to separate multiples ones of the co-channel beacon transmissions in the snapshot or eliminate inter-symbol interference and inter-subcarrier interference in the snapshot.

2. The method of claim 1, wherein at least one of generating or using is configured to be performed at a network user, a network node, or a network operations center.

3. The method of claim 2, wherein the network operations center employs a wireless communication link configured to provision each of a plurality of network nodes with configuration data or time symbols.

4. The method of claim 2, wherein generating is performed at the network user or the network node, the method further comprising communicating the snapshot to the network operations center, the network operations center being configured to perform the subcarrier demodulation, code nulling, or Class-C linear MMSE operation.

5. The method of claim 3, wherein generating and using are performed at the network user, the method further comprising: configuring the network user to receive network node locations and configuration data from the network operations center via the wireless communication link; and configuring the network user to process the snapshot with the network node locations and configuration data for determining position and timing.

6. The method of claim 2, wherein generating is performed at the network user, and wherein the snapshot contains co-channel beacon transmissions transmitted from a plurality of network nodes within the network user's field of view.

7. The method of claim 3, wherein the wireless communication link comprises Long Term Evolution (LTE), fourth-generation wireless (4G), fifth generation (5G) new radio (NR), IEEE 802.11 wireless local area network (WLAN), Zigbee, or Bluetooth technology.

8. The method of claim 1, wherein the snapshot comprises a prefix and a suffix to encompass inter-slot interference introduced by timing offset between the plurality of beacon transmitters and the receiver.

9. The method of claim 1, further comprising channelizing the superposition to separate the co-channel beacon transmissions in the snapshot into frequency subcarriers and time symbols to produce a channelized snapshot; wherein using is configured to perform the subcarrier demodulation, code nulling, or Class-C linear MMSE operation on the channelized snapshot.

10. The method of claim 9, further comprising stacking the channelized snapshot into an N.sub.data×M.sub.DoF windowed data matrix, wherein M.sub.DoF denotes degrees-of-freedom (DoF) of the windowed data matrix, each of the plurality of co-channel beacon transmissions being redundant across the DoF, and wherein N.sub.data denotes a number of data samples in the windowed data matrix.

11. The method of claim 1, wherein using is configured to separate the co-channel beacon transmissions with precision dictated by received power of each of the co-channel beacon transmissions above a receiver noise floor, and irrespective of other ones of the co-channel beacon transmissions received at the same time and frequency.

12. The method of claim 1, further comprising determining geo-observables from separated co-channel beacon transmissions.

13. The method of claim 12, further comprising: determining at least one of positioning and timing from the geo-observables.

14. An apparatus, comprising: a receiver configured to receive a superposition of co-channel beacon transmissions transmitted from a plurality of beacon transmitters, each of the co-channel beacon transmissions having at least one of spectral redundancy and temporal redundancy; the receiver configured for generating a snapshot of the superposition; and a position/timing (P/T) solution generator communicatively coupled to the receiver, the P/T solution generator configured to employ subcarrier demodulation, code nulling, or a Class-C linear minimum-mean-square error (MMSE) operation to separate multiples ones of the co-channel beacon transmissions in the snapshot or eliminate inter-symbol interference and inter-subcarrier interference in the snapshot.

15. The apparatus of claim 14, wherein the P/T solution generator is located in a network user, a network node, or a network operations center.

16. The apparatus of claim 15, wherein the network operations center employs a wireless communication link configured to provision each of a plurality of network nodes with configuration data or time symbols.

17. The apparatus of claim 15, wherein the receiver is located at the network user or the network node, and the P/T solution generator is located at the network operations center, the apparatus further comprising a communication network configured for communicating the snapshot to the network operations center.

18. The apparatus of claim 16, wherein the receiver and the P/T solution generator are located at the network user, the receiver being configured to receive network node locations and configuration data from the network operations center via the wireless communication link; and the P/T solution generator being configured to process the snapshot with the network node locations and configuration data for determining position and timing.

19. The apparatus of claim 15, wherein the receiver is located at the network user, and wherein the snapshot contains co-channel beacon transmissions transmitted from a plurality of network nodes within the network user's field of view.

20. The apparatus of claim 16, wherein the wireless communication link comprises Long Term Evolution (LTE), fourth-generation wireless (4G), fifth generation (5G) new radio (NR), IEEE 802.11 wireless local area network (WLAN), Zigbee, or Bluetooth technology.

21. The apparatus of claim 14, wherein the snapshot comprises a prefix and a suffix to encompass inter-slot interference introduced by timing offset between the plurality of beacon transmitters and the receiver.

22. The apparatus of claim 14, wherein the receiver comprises a channelizer configured for channelizing the superposition to separate the co-channel beacon transmissions in the snapshot into frequency subcarriers and time symbols to produce a channelized snapshot; wherein the P/T solution generator is configured to perform the subcarrier demodulation, code nulling, or Class-C linear MMSE operation on the channelized snapshot.

23. The apparatus of claim 22, wherein the P/T solution generator is configured for stacking the channelized snapshot into an N.sub.data×M.sub.DoF windowed data matrix, wherein M.sub.DoF denotes degrees-of-freedom (DoF) of the windowed data matrix, each of the plurality of co-channel beacon transmissions being redundant across the DoF, and wherein N.sub.data denotes a number of data samples in the windowed data matrix.

24. The apparatus of claim 14, wherein the P/T solution generator is configured to separate the co-channel beacon transmissions with precision dictated by received power of each of the co-channel beacon transmissions above a receiver noise floor, and irrespective of other ones of the co-channel beacon transmissions received at the same time and frequency.

25. The apparatus of claim 14, wherein the P/T solution generator is configured for determining geo-observables from separated co-channel beacon transmissions.

26. The apparatus of claim 25, wherein the P/T solution generator is configured for determining at least one of positioning and timing from the geo-observables.

27. An apparatus configured to receive a superposition of co-channel beacon transmissions transmitted from a plurality of beacon transmitters, each of the co-channel beacon transmissions having at least one of spectral redundancy and temporal redundancy; the apparatus comprising: at least one processor and at least one memory in electronic communication with the at least one processor, the at least one memory having instructions stored therein and executable by the at least one processor for: generating a snapshot of the superposition; and using subcarrier demodulation, code nulling, or a Class-C linear minimum-mean-square error (MMSE) operation to separate multiples ones of the co-channel beacon transmissions in the snapshot or eliminate inter-symbol interference and inter-subcarrier interference in the snapshot.

28. The apparatus of claim 27, wherein the at least one processor and the at least one memory are located at a network user, a network node, or a network operations center.

29. The apparatus of claim 28, wherein the network operations center employs a wireless communication link configured to provision each of a plurality of network nodes with configuration data or time symbols.

30. The apparatus of claim 28, wherein generating is performed at the network user or the network node, the apparatus further comprising a transmitter configured for communicating the snapshot to the network operations center, the network operations center being configured to perform the subcarrier demodulation, code nulling, or Class-C linear MMSE operation.

31. The apparatus of claim 29, wherein generating and using are performed at the network user, and the at least one memory has instructions stored therein and executable by the at least one processor for: receiving network node locations and configuration data from the network operations center via the wireless communication link; and processing the snapshot with the network node locations and configuration data for determining position and timing.

32. The apparatus of claim 28, wherein generating is performed at the network user, and wherein the snapshot contains co-channel beacon transmissions transmitted from a plurality of network nodes within the network user's field of view.

33. The apparatus of claim 29, wherein the wireless communication link comprises Long Term Evolution (LTE), fourth-generation wireless (4G), fifth generation (5G) new radio (NR), IEEE 802.11 wireless local area network (WLAN), Zigbee, or Bluetooth technology.

34. The apparatus of claim 27, wherein the snapshot comprises a prefix and a suffix to encompass inter-slot interference introduced by timing offset between the plurality of beacon transmitters and the receiver.

35. The apparatus of claim 27, wherein the at least one memory has instructions stored therein and executable by the at least one processor for: channelizing the superposition to separate the co-channel beacon transmissions in the snapshot into frequency subcarriers and time symbols to produce a channelized snapshot; wherein using is configured to perform the subcarrier demodulation, code nulling, or Class-C linear MMSE operation on the channelized snapshot.

36. The apparatus of claim 35, wherein the at least one memory has instructions stored therein and executable by the at least one processor for: stacking the channelized snapshot into an N.sub.data×M.sub.DoF windowed data matrix, wherein M.sub.DoF denotes degrees-of-freedom (DoF) of the windowed data matrix, each of the plurality of co-channel beacon transmissions being redundant across the DoF, and wherein N.sub.data denotes a number of data samples in the windowed data matrix.

37. The apparatus of claim 27, wherein using is configured to separate the co-channel beacon transmissions with precision dictated by received power of each of the co-channel beacon transmissions above a receiver noise floor, and irrespective of other ones of the co-channel beacon transmissions received at the same time and frequency.

38. The apparatus of claim 27, wherein the at least one memory has instructions stored therein and executable by the at least one processor for: determining geo-observables from separated co-channel beacon transmissions.

39. The apparatus of claim 38, wherein the at least one memory has instructions stored therein and executable by the at least one processor for: determining at least one of positioning and timing from the geo-observables.

40. A non-transitory computer-readable memory, and instructions stored therein and executable by at least one processor for: receiving a superposition of co-channel beacon transmissions transmitted from a plurality of beacon transmitters, each of the co-channel beacon transmissions having at least one of spectral redundancy and temporal redundancy; generating a snapshot of the superposition; and using subcarrier demodulation, code nulling, or a Class-C linear minimum-mean-square error (MMSE) operation to separate multiples ones of the co-channel beacon transmissions in the snapshot or eliminate inter-symbol interference and inter-subcarrier interference in the snapshot.

41. The non-transitory computer-readable memory of claim 40, wherein at least one of receiving, generating, or using are configured to be performed at a network user, a network node, or a network operations center.

42. The non-transitory computer-readable memory of claim 40, further comprising instructions stored therein and executable by the at least one processor for employing a wireless communication link to provision each of a plurality of network nodes with configuration data or time symbols.

43. The non-transitory computer-readable memory of claim 41, further comprising instructions stored therein and executable by the at least one processor for communicating the snapshot to the network operations center, the network operations center being configured to perform the subcarrier demodulation, code nulling, or Class-C linear MMSE operation.

44. The non-transitory computer-readable memory of claim 42, further comprising instructions stored therein and executable by the at least one processor for: receiving network node locations and configuration data from the network operations center via the wireless communication link; and processing the snapshot with the network node locations and configuration data for determining position and timing.

45. The non-transitory computer-readable memory of claim 41, wherein generating is performed at the network user, and wherein the snapshot contains co-channel beacon transmissions transmitted from a plurality of network nodes within the network user's field of view.

46. The non-transitory computer-readable memory of claim 42, wherein the wireless communication link comprises Long Term Evolution (LTE), fourth-generation wireless (4G), fifth generation (5G) new radio (NR), IEEE 802.11 wireless local area network (WLAN), Zigbee, or Bluetooth technology.

47. The non-transitory computer-readable memory of claim 40, wherein the snapshot comprises a prefix and a suffix to encompass inter-slot interference introduced by timing offset between the plurality of beacon transmitters and the receiver.

48. The non-transitory computer-readable memory of claim 40, further comprising instructions stored therein and executable by the at least one processor for: channelizing the superposition to separate the co-channel beacon transmissions in the snapshot into frequency subcarriers and time symbols to produce a channelized snapshot; wherein using is configured to perform the subcarrier demodulation, code nulling, or Class-C linear MMSE operation on the channelized snapshot.

49. The non-transitory computer-readable memory of claim 48, further comprising instructions stored therein and executable by the at least one processor for: stacking the channelized snapshot into an N.sub.data×M.sub.DoF windowed data matrix, wherein M.sub.DoF denotes degrees-of-freedom (DoF) of the windowed data matrix, each of the plurality of co-channel beacon transmissions being redundant across the DoF, and wherein N.sub.data denotes a number of data samples in the windowed data matrix.

50. The non-transitory computer-readable memory of claim 40, wherein using is configured to separate the co-channel beacon transmissions with precision dictated by received power of each of the co-channel beacon transmissions above a receiver noise floor, and irrespective of other ones of the co-channel beacon transmissions received at the same time and frequency.

51. The non-transitory computer-readable memory of claim 40, further comprising instructions stored therein and executable by the at least one processor for: determining geo-observables from separated co-channel beacon transmissions.

52. The non-transitory computer-readable memory of claim 51, further comprising instructions stored therein and executable by the at least one processor for: determining at least one of positioning and timing from the geo-observables.

Description

BRIEF DESCRIPTION OF THE DRAWINGS

[0061] A further understanding of the nature and advantages of the present disclosure may be realized by reference to the following drawings. In the appended figures, similar components or features may have the same reference label.

[0062] FIG. 1 illustrates a resilient distributed transmitter network that can be instantiated according to some aspects of the disclosure.

[0063] FIG. 2 illustrates an alternate resilient distributed receiver network that can be instantiated according to some aspects of the disclosure.

[0064] FIG. 3 illustrates an alternate resilient distributed transponder network that can be instantiated according to some aspects of the disclosure.

[0065] FIGS. 4A and 4B depict exemplary network scenarios in which network users are 100 Class-1 sUAS's.

[0066] FIGS. 5A and 5B shows the received incident power (RIP) of each link in the network scenarios shown in FIG. 4A and FIG. 4B, as a function of the time-of-arrival (TOA) and frequency-of-arrival (FOA) of beacon transmitters in the field of view (FoV) of each receiver.

[0067] FIG. 6A and FIG. 6B shows the observed SINR of each link in the network scenarios shown in FIG. 4, as a function of the beacon TOA and FOA observed at a representative beacon receiver.

[0068] FIG. 7 and FIG. 8 depicts operations used to generate and transmit a beacon.

[0069] FIG. 9A depicts spectral redundancy in a beacon when a Kronecker product operation is used to construct a full subcarrier vector from simple BPSK inner subcarrier vectors and BPSK outer subcarrier vectors; and FIG. 9B illustrates temporal structure of the beacon for an aspect in which time symbols are organized into slots.

[0070] FIG. 10A-10D show the frequency and time response of two beacons generated in accordance with an aspect of the disclosure.

[0071] FIG. 11 shows a user receiver system concept used in an aspect of the disclosure.

[0072] FIG. 12 shows exemplary reception and PNT processing operations used in “cold start” scenarios where the receiver has no prior knowledge of its position or precise timing and carrier synchronization to the beacon transmissions.

[0073] FIG. 13 depicts data channelization operations employed by aspects of the disclosure.

[0074] FIG. 14 shows the frequency response of an exemplary channelizer window developed in accordance with some aspects of the disclosure.

[0075] FIG. 15 illustrates an efficient resilient TOA-FOA surface generation process.

[0076] FIG. 16A and FIG. 16B and FIG. 17A and FIG. 17B show the TOA-FOA surface obtained for two aspects of the disclosure.

[0077] FIG. 18, FIG. 19, and FIG. 20 depict geo-observable based positioning/timing operations used in an aspect of the disclosure.

[0078] FIG. 21A-FIG. 21C illustrate joint geolocation and synchronization procedure for an exemplary reception scenario.

[0079] FIG. 22A-FIG. 22C and FIG. 23A-FIG. 23C show performance of the method for all of the network users in two exemplary reception scenarios.

DETAILED DESCRIPTION OF THE DRAWINGS

[0080] The detailed description set forth below in connection with the appended drawings is intended as a description of various configurations and is not intended to represent the configurations in which the concepts described herein may be practiced. The detailed description includes specific details for the purposes of providing a thorough understanding of various concepts. However, it will be apparent to those skilled in the art that these concepts may be practiced without these specific details. In some instances, well-known structures and components are shown in block diagram form in order to avoid obscuring such concepts.

Exemplary Resilient Distributed Positioning Networks

[0081] FIG. 1 illustrates a resilient distributed transmitter network (RDTN) that can be instantiated using the system introduced here. In this network, beacons are transmitted from a network of L network nodes 401 on a known frequency channel. Exemplary network nodes 401 include fixed outdoor transmitters co-located with cellular transmission towers or 802.11 access points; indoor transmitters coexisting with 802.11 WLAN or 802.15 Bluetooth or Zigbee networks; or standalone transmitters. The network nodes 401 are connected to a network operations center (NOC) 403, e.g., a 5GNR multi-access edge computer (MEC), or a UAS service supplier (USS), which provisions each of those network nodes 401 with separate configuration data or time symbols over a secure communication link. In alternate instantiations, the network nodes 401 choose their own configuration and apprise the NOC 403 site of that choice. The transmitted beacons are received by M users 400, e.g., Class-1 sUAS's, as shown in this FIG., ground vehicles, cellular user equipment (UE's), 802.11 STA's, 802.15 Zigbee or Bluetooth devices, or Internet of Thing devices, etc., each of which collect snapshots of data containing a superposition of the beacons transmit from network nodes 401 within their field of view (FoV). In one aspect, the users 400 backhaul those snapshots to the NOC 403 over at least one wireless communication transceiver 402, e.g., a cellular LTE, 4G, or 5G network or an 802.11 wireless local-area network (WLAN), connected to the users 400 and the NOC 403. In alternate aspects, the users 400 may perform positioning and timing operations themselves, using network node 401 locations and configuration data provided through the wireless communication transceiver 402, and transmit results of those operations to the NOC 403 over the wireless communication transceiver 402. In some aspects, different users 400 with varying capabilities and network permissions, including time-varying capabilities and permissions, may implement one aspect or the other.

[0082] FIG. 2 illustrates an alternate resilient distributed receiver network (RDRN) that can be instantiated using the system introduced here. In this network, beacons are transmitted from the users 400, based on beacon configuration data transmitted to the users 400 from the NOC 403 over a wireless communication transceiver 402. Snapshots can be transmitted on a continuous or user-400 initiated basis, e.g., after detection of PNT outage conditions, or upon request by the NOC 403. The snapshots are then received at network nodes 401, and backhauled to the NOC 403 from the network nodes 401.

[0083] FIG. 3 illustrates an alternate resilient distributed transceiver network (RDXN) that can also be instantiated using the system introduced here. In this network, beacons are transmitted from network users 401, using beacon time symbols or configuration data provided by NOC 403 and are received and directly transponded from the users 400 as needed by the users 400 or the NOC 403.

Exemplary Class-1 sUAS Deployment Scenarios

[0084] FIG. 4A and FIG. 4B depicts two exemplary resilient distributed positioning network deployment scenarios, one intended for operation in Subband A of the Location and Monitoring Services (LMS) Band, located at 904-906 MHz (FIG. 4A), referred to here as the “LMS Scenario”, and the other intended for operation in 802.11 2.4 GHz Channel 13, located at 2,462-2,482 MHz (FIG. 4B), a 20 MHz channel in the 2.4-2.485 ISM band that is not used for 802.11 transmissions in the United States, referred to here as the “2.4 GHz Ch. 13 Scenario”. In each scenario, 100 users 400 are deployed at 3-to-122 meter (10-400 foot) elevations above ground level (AGL), and 22.3-to-44.7 meter/second (50-100 mile/hour) airspeeds, consistent with class-1 sUAS's, over a 32.2 km (20 mile) range from an areal analysis center. Assuming spherical Earth model with 4/3 refraction index, the users 400 have a 45.5 km (28.3 mile) ground horizon at maximum altitude. In FIG. 4A the users 400 are in the FoV of 32 fixed network nodes 401, deployed at 15.2-to-45.7 (50-150 foot) elevations AGL, i.e., with a ground horizon of 27.9 km (17.3 miles) at maximum altitude, and in a roughly hexagonal layout, consistent with sparse deployment from cell towers, over a 17,227 square-kilometer area. In FIG. 4B the users 400 are in the FoV of 50 fixed network nodes 401, deployed over the same range of altitudes, and over a 16,480 square-kilometer area. Both layouts cover the full FoV of the users 400 over the range of users 400 and network node 401 elevations, allowing 1,382 links in the LMS Scenario, or between 4 and 25 links per user 400; and 2,517 links in the 2.4 GHz Ch. 13 Scenario, or between 5 and 41 links per user 400.

[0085] FIG. 5A and FIG. 5B show the received incident power (RIP) of the beacons at each receiver, as a function of the time-of-arrival (TOA) and frequency-of-arrival (FOA), i.e., Doppler shift, of those received beacons, for the LMS Scenario (FIG. 5A) and the 2.4 GHz Ch. 13 Scenario (FIG. 5B). In each scenario, the beacons are assumed to possess azimuthally omnidirectional antennas with 6 dBi gain at 0° elevation and −3 dBi gain beyond ±300 elevation. In the LMS Scenario, the beacons are transmitted at a 40.9 dBm power, yielding an ERP of 30 Watts, the maximum power allowed in this band. In the 2.4 GHz Ch. 13 Scenario, the beacons are transmitted at 30 dBm power, yielding a 36 dBm EIRP, the maximum FCC Part 15 compliant power allowed in the 2.4-2.485 GHz ISM band. Each link assumes a two-ray median pathloss, with 4 dB rms log-normal shadowing and (20 dB mean, 4 dB standard deviation) Rician fading over each link.

[0086] As these FIGS. show, the link TOA's are restricted to 3-229 μs in the LMS Scenario, and between 4-232 μs in the 2.4 GHz Ch. 13 Scenario, much lower than the 67-94 ms TOA range expected for GNSS signals transmitted from MEO. Similarly, the link FOA's range are restricted to ±135 Hz in the LMS Scenario and ±369 Hz in the 2.4 GHz Ch. 13 Scenario, much lower than the ±6 kHz FOA range expected for L-band GNSS signals transmitted from MEO. This is an exploitable feature of the RDPN for both scenarios. At the same, the RIP of the beacons range from −88 dBm to −40 dBm in the LMS Scenario, and from −105 dBm to −64 dBm in the 2.4 GHz Ch. 13 Scenario, much stronger than the nominal −130 dBm GPS L1 signal strength at the Earth surface. While this is also a clear advantage for any beacon-based positioning solution, it also shows that the beacons will likely be received at positive SNR and with significant near-far interference. That is, the performance of conventional “matched filter” receivers that correlate the received signal against replicas of the transmitted beacons will be limited by the interference observed relative to each beacon, despite their high receive SNR, due to self-interference between those received beacons.

[0087] This observation is borne out in FIG. 6A and FIG. 6B, which show the observed TOA, FOA, and SINR of each link shown in FIG. 5A and FIG. 5B, for the LMS Scenario (FIG. 6A) and the 2.4 GHz Ch. 13 Scenario (FIG. 6B). In both cases, the beacon transmitters assume azimuthally omnidirectional antennas with 6 dBi gain at 0° elevation and −3 dBi gain beyond ±300 elevation, and receivers with 6 dB rolled-up noise figure (NF) through the receiver's analog-to-digital converter (ADC), and assume the ADC's and downconversion local oscillators (LO) are locked to clocks with ±10 ppm and ±100 μs 90th percentile rate and timing offset, respectively. In this scenario, the beacon transmitters are assumed to have negligible timing and LO offset, i.e., they are employing GPS disciplined oscillators, or they are synchronized to UTC via a calibration procedure.

[0088] As these FIGS. show, the ADC output SINR ranges between −46 dB and +25 dB in the LMS Scenario, and between −39 dB and +8 dB in the 2.4 GHz Ch. 13 Scenario. In fact, while all of the links are above 0 dB SNR in the LMS Scenario, only 3.5% of the links are above a 0 dB SINR, and less than 18% of the links are above a −10 dB SINR. Similarly, 82% of the links are above 0 dB SNR in the 2.4 GHz Ch. 13 Scenario, but only 1.2% of the links are above 0 dB SINR, and 10% of the links are above a −10 dB SINR. Hence, the received beacons are clearly in an interference-limited environment. This is the reason that competing systems introduce time hopping and time slotting into their beacon transmitters—in order to avoid such interference.

[0089] These results motivate the use of beacons that can both exploit the much tighter range of TOA and FOA obtaining in a ground-to-low-altitude reception geometry, and allow the use of interference excision methods that can separate the received beacons with performance gated by their (high) receive SNR, rather than their (low) receive SINR.

Beacon Generation and Transmission System

[0090] FIG. 7 depicts operations used to generate beacon time symbols for the networks shown in FIG. 1 and FIG. 3, and for aspects in which the network nodes 401 are provided beacon time symbols from the NOC 403. The time symbol generator first selects code indices 101 for each network node 401, creating inner and outer code indices 101 {c.sub.0(),c.sub.1, where the code indices 101 point to code libraries 102 containing inner and outer code vectors {θ(c.sub.0)}.sub.c.sub.0.sub.=0.sup.c.sup.0.sup.−1 and {θ(c.sub.1)}.sub.c.sub.1.sub.=0.sup.c.sup.1.sup.−1, respectively, and where θ(c.sub.i)=[θ.sub.k.sub.i(c.sub.i)].sub.k.sub.i.sub.=0.sup.K.sup.i.sup.−1 is a K.sub.i×1 vector of phase values in radians. The code indices 101 can be chosen deterministically, pseudorandomly, or unpredictably, or based on environmental factors determined by the network, e.g., evidence of spoofing by malicious actors.

[0091] In one aspect, the phase vectors contained in the code libraries 102 are designed to yield library codebook signal

[00001]$s\left(t,\theta \left(c_i\right)\right)=\frac{1}{\sqrt{K_i}}\underset{k_i=0}{\overset{K_i-1}{.Math.}}\exp \left(j\left({\theta }_{k_i}\left(c_i\right)+2\pi k_{i}t\right)\right)$

with low peak-to-average power ratio (PAPR), and low cross-correlation between other library codebook signals. The code libraries 102 can be designed to satisfy other system requirements as well.

[0092] The time symbol generator then performs an inner subcarrier construction operation 103 to generate K.sub.0×1 inner subcarrier vector 104 =[(k.sub.0)].sub.k.sub.0.sub.=0.sup.K.sup.0.sup.−1 and an outer subcarrier construction operation 105 to generate K.sub.1×1 outer subcarrier vector 106 =[(k.sub.1)].sub.k.sub.1.sub.=0.sup.K.sup.1.sup.−1. In one aspect, the vectors are generated using simple phase modulation operations, such that

[00002]${S_{T_0\left(\ell \right)}=\frac{1}{\sqrt{K_0}}[\exp (j{\theta }_{T_0\left(\ell \right)}\left(k_0\right)\}]}_{k_0=0}^{K_0-1}\mathrm{and}{S}_{T_1\left(\ell \right)}={\frac{1}{\sqrt{K_1}}[\exp (j{\theta }_{T_1\left(\ell \right)}\left(k_1\right)\}]}_{k_1=0}^{K_1-1}\mathrm{where}{\theta }_{T_i\left(\ell \right)}\left(k_i\right)={\theta }_{k_i}\left(c_i\left(\ell \right)\right).$

Other aspects apply non-uniform amplitude weightings to either or both vectors, for example, to further reduce PAPR of the transmitted beacon, reduce interference to non-beacon networks caused by beacons in selected portions of the beacon transmission band, or reduce susceptibility to non-beacon interference at network receivers operating in the beacon transmission band.

[0093] The inner subcarrier vector 104 and outer subcarrier vector 106 are then combined 107 to form K.sub.sub×1 full subcarrier vector 108 =[(k.sub.sub)].sub.k.sub.sub.sub.=0.sup.K.sup.sub.sup.−1. In one aspect, this is performed using a Kronecker product operation, such that =.Math.. In other aspects, this may be a more complex combining operation, for example, to improve robustness to LO frequency uncertainty at the transmitter or receiver, reduce interference to non-beacon networks caused by beacons in selected portions of the beacon transmission band, or reduce susceptibility to non-beacon interference at network receivers operating in the beacon transmission band.

[0094] The full subcarrier vector 108 is then passed through an optional subcarrier preemphasis 106 operation to generate preemphasized subcarrier vector =, where “∘” denotes the element-wise (Shur or Hadamard) product operation, and is a K.sub.sub×1 preemphasis vector 110 that compensates for front-end digital-to-analog conversion (DAC), antialiasing lowpass filtering (LPF), and upconversion operations performed at network node 401 . The preemphasis vector 110 can be designed using analytic models for beacon transmission operations 213; or using calibration data obtained at each network node 401, for example, as described in G. Pattabiraman, S. Melyappan, A. Raghupathy, H. Sankar, “Wide Area Positioning System,” U.S. Pat. No. 8,130,141, issued March 2012, and can be based on the magnitude or complex value of those beacon transmission operations 213.

[0095] The preemphasized subcarriers are then passed to a multitone modulator 111 that transforms the subcarriers to the time domain, and is optionally quantized 112, to provide an N.sub.DAC×1 time symbol vector 113 =[(n.sub.DAC)].sub.n.sub.DAC.sub.=0.sup.N.sup.DAC.sup.−1. The multitone modulator 111 can be implemented using combinations of discrete Fourier transform (DFT), inverse DFT (IDFT), fast Fourier transform (FFT), or inverse FFT (IFFT) operations, or using polyphase filtering or discrete filter bank operations.

[0096] Each time symbol vector 113 is then passed from the NOC 403 to the beacon transmitter over a beacon communication bus 114. Exemplary communication networks supporting a beacon communication bus 114 can include Ethernet-based networks, optical networks, power-line communication (PLC) networks, 802.11 WLAN's, 802.15 Zigbee or Bluetooth networks, or 3G, 4G LTE, or 5G cellular networks. In the networks shown in FIG. 1 and FIG. 3, including this aspect of the disclosure, the beacon communication bus 114 connects the NOC 403 to the network nodes 401. In the network shown in FIG. 2, the beacon communication bus 114 connects the NOC 403 to the users 400 over a wireless communication transceiver 402.

[0097] FIG. 8 shows operations performed in this aspect of the disclosure at each beacon transmitter, e.g., the network nodes 401 for the networks shown in FIG. 1 and FIG. 3, or the users 400 for the network shown in FIG. 2. At each beacon transmitter, the N.sub.DAC×1 time symbol vector 113 is first obtained from the beacon communication bus 114. If K.sub.sub is even, the time symbol vector 113 is placed through time symbol extension operation 211 to convert the N.sub.DAC×1 time symbol vector 113 to a 2N.sub.DAC×1 time symbol vector 113 given by

[00003]${\tilde{s}}_{T\left(\ell \right)}\leftarrow \left(\begin{array}{c}{\tilde{s}}_{T\left(\ell \right)}\\ -{\tilde{s}}_{T\left(\ell \right)}\end{array}\right),$

and placed in local storage 212. For the networks shown in FIG. 1 and FIG. 3, the time symbol vector 113 in local storage 212 is then continuously and repeatedly cycled through beacon transmission operations 213, until a new time symbol is provided by the NOC 403, or it must cease transmission for other reasons, e.g., to comply with license provisions of spectra reserved for time-division-duplex (TDD) operation, or based on other instructions received from the NOC 403. For the network shown in FIG. 2, the time symbol vector 113 is transmitted from the user 400 over varying time interval durations, dependent for example upon the priority of positioning/timing solutions demanded by the user 400 or the NOC 403, or user 400 energy conservation needs.

[0098] The operations used to perform the beacon transmission operations 213 are analogous to an arbitrary waveform generator (AWG). In the aspect shown in FIG. 8, the operations comprise dual DAC, dual LPF, frequency shift operations using in-phase and quadrature LO's, and power amplification (PA) operations. However, many other means can be used to implement the beacon transmission operations 213, e.g., superheterodyne transmitters that convert the time symbol vector 113 to a real intermediate-frequency (real-IF) form and perform the frequency conversion in multiple steps, and polar modulators, for example, D. Kirkpatrick, E. McCune, Jr, “Polar Modulation Using Product Mode,” U.S. Pat. No. 9,397,713, issued July 2016, which convert the time symbol vector 113 to polar form and separately phase modulate the time symbol vector 113 phase and amplitude-modulate the time symbol vector 113 amplitude. These conversions can be implemented at the NOC 403, or at the beacon transmitters.

[0099] Typically, the DAC and LO employed in the beacon transmission operations 213 are locked to a system clock 214, which in general, has a clock rate and timing that is offset from a common time standard, e.g., UTC. As shown in FIG. 8, consistent with the networks shown in FIG. 1 and FIG. 3, the system clock 214 for network node 401 is offset from UTC by rate offset and timing offset =−t.sub.UTC(), where is the internal clock time at UTC time t.sub.UTC(). In some aspects, the system clock 214 is synchronized to an external time and frequency standard using an external source, e.g., a GNSS receiver (Rx) 215. In some aspects, the system clock 214 is brought into a common time standard using network calibration methods computed at the NOC 403, e.g., by providing the system clock 214 with clock synchronization data 215, e.g., timing and rate offset estimates and . In the network shown in FIG. 2, the system clock 214 employed each user 400 can similarly be brought into synchronization with a common time standard using clock synchronization data 215 provided by the NOC 403.

[0100] In some aspects consistent with the networks shown in FIG. 1 and FIG. 3, the network nodes 401 are provisioned with the code library 101, and the NOC 403 sends it the code indices 101 to be used to generate time symbol vector 113. Alternately, each network node 401 can select its own code indices 101 and communicate it back to the NOC 403. These aspects greatly reduce the data rate needed to communicate between the network nodes 401 and NOC 403, but increase operation security of the network, by requiring the beacons to possess copies of the code libraries 102.

[0101] For an exemplary subcarrier frequency layout assumed here, FIG. 9A illustrates the spectral redundancy imposed in aspects where a Kronecker product operation is used to construct the full subcarrier vector 108 from simple BPSK inner subcarrier vectors 102 and BPSK outer subcarrier vectors 103; and FIG. 9B illustrates the temporal structure of the beacon for an aspect in which the time symbols are organized into slots with N.sub.rep per slot. In FIG. 9A, the inner subcarrier vector 103 is a 5-element vector [+1 +1 −1 −1 +1].sup.T, while the outer subcarrier vector 106 is a 4-element vector [+1 −1 −1 +1].sup.T. The inner and outer subcarrier combining operation 107 then yields a 20-element full subcarrier vector 108 in which the 5-phase inner subcarrier vector 103 is replicated over 4 clusters, each of which is modulated by one element of the outer subcarrier vector 106. Beacon generated by the multitone modulator 107 then has complex baseband representation

[00004]$\begin{array}{cc}\begin{array}{c}{s}_{T\left(\ell \right)}\left(t\right)=\underset{k_{sub}=0}{\stackrel{K_{sub}-1}{.Math.}}{S}_{T\left(\ell \right)}\left(k_{sub}\right)\exp \left\{j2\pi f\left(k_{sub}\right)t\right\}\\ =\frac{1}{\sqrt{K_{sub}}}\underset{k_{sub}=0}{\overset{K_{sub}-1}{.Math.}}\exp (j\left({\theta }_{T\left(\ell \right)}\left(k_{sub}\right)+2\pi f\left(k_{sub}\right)t\right)\}\square \end{array}& \left(\mathrm{Eq1}\right)\end{array}$

at the PA input in the beacon transmission operations 213, where the subcarrier frequencies are given by

[00005]$\begin{array}{cc}f\left(k_{sub}\right)=\left(k_{sub}-\frac{K_{sub}-1}{2}\right)f_{sym},k_{sub}=0,\rightleftharpoons ,K_{sub}-1\square & \left(\mathrm{Eq2}\right)\end{array}$

and where f.sub.sym=1/T.sub.sym is the subcarrier spacing. Assuming the internal subcarrier structure shown here, then (t)=(t), where

[00006]$\begin{array}{cc}{s}_{T_i\left(\ell \right)}\left(t\right)=\frac{1}{\sqrt{K_i}}\underset{k_i=0}{\overset{K_i-1}{.Math.}}\exp (j\left({\theta }_{T_i\left(\ell \right)}\left(k_i\right)+2\pi f_i\left(k_i\right)t\right)\}\square & \left(\mathrm{Eq3}\right)\end{array}$

and where inner and outer subcarrier frequencies f.sub.0(k.sub.0) and f.sub.1(k.sub.1) are given by

[00007]$\begin{array}{cc}{f}_0\left(k_0\right)=\left(k_0-\frac{K_0-1}{2}\right)f_{sym},k_0=0,\rightleftharpoons ,K_0-1,& \left(\mathrm{Eq4}\right)\end{array}$$\begin{array}{cc}{f}_1\left(k_1\right)=\left(k_1-\frac{K_1-1}{2}\right)K_0{f}_{\mathrm{sym}},k_1=0,\rightleftharpoons ,K_1-1\square & \left(\mathrm{Eq5}\right)\end{array}$

respectively, such that K.sub.sub=K.sub.0K.sub.1, f(K.sub.0k.sub.1+k.sub.0)=f.sub.0(k.sub.0)+f.sub.1(k.sub.1) and (K.sub.0k.sub.1+k.sub.0)=(k.sub.0)=(k.sub.1).

[0102] This beacon can be interpreted as a stacked-carrier spread spectrum (SCSS) signal, in which a narrowband signal (t) with bandwidth K.sub.0f.sub.sym and period 2.sup.(K.sup.0.sup.−1) mod 2T.sub.sym is spread by a wideband signal (t) with bandwidth K.sub.subf.sub.sym and period 2.sup.(K.sup.1.sup.−1) mod 2T.sub.sym/K.sub.0. The resultant beacon possesses massive spectral redundancy, both between clusters (replication of inner code over K.sub.1 independent clusters), and within clusters (replication of outer code over K.sub.0 subcarriers within each cluster).

[0103] FIG. 9B further illustrates the additional temporal redundancy in the beacon over each symbol repetition. If the time symbol possesses an even number of subcarriers, then the subcarriers will be offset from the frequency origin by a factor of f.sub.sym/2, similar to LTE SC-FDMA uplink signals. In this case, successive time symbol repetitions needed to be inverted to preserve signal energy within each subcarrier. If the time symbol possesses an odd number of subcarriers, then this offset is removed, similar to LTE OFDM downlink signals, and this successive inversion is not needed. The replication induces temporal redundancy in either case.

Exemplary Class-1 sUAS Compatible Beacon Generation and Transmission Parameters

[0104] Table 1 lists exemplary beacon generation and transmission parameters compatible with TOA and FOA ranges expected for Class-1 sUAS's, and for network geometries shown in FIG. 4. The key parameter providing compatibility with Class-1 sUAS's is the 250 μs symbol duration T.sub.sym, which encompasses the full range of TOA's shown in FIG. 5, driven by the maximum altitude allowed for Class-1 sUAS's. The resultant 4 kHz subcarrier frequency separation is also compatible with the range of FOA's shown in FIG. 5, driven by the maximum airspeed allowed for Class-1 sUAS's. In each case, the beacon bandwidth fits into the channel allocation for each band. The 80-element inner subcarrier vector 104 dimension chosen for the 2.4 GHz Ch. 13 Scenario allows up to 80 co-channel beacons to be separated using linear-algebraic signal separation methods. The 24-element inner subcarrier vector 104 dimension chosen for the LMS Scenario—driven largely by the much narrower 2 MHz bandwidth constraint in this subband—only allows up to 24 co-channel beacons to be separated using linear-algebraic methods, hence a sparser set of co-channel beacons can be deployed. The narrower bandwidth and sparser deployment degrades the geolocation capability of this network. However, the higher transmit power requirement compensates for much of this performance loss. It should be noted that correlative techniques would be strongly interference-limited in this case, and could not exploit any of the power advantage available in this subband.

TABLE-US-00001 TABLE 1 Exemplary Class-1 sUAS Compatible Beacon Transmission Parameters Parameter LMS Scenario 2.4 GHz Ch. 13 Scenario f.sub.T 905 MHz 2,472 MHz P.sub.T 40.9 dBm 30 dBm G.sub.T, 0° elevtion 6 dBi (30 W ERP) 6 dbi (36 dBm EIRP) T.sub.sym 250 μs (4 kHz f.sub.sym) 250 μs (4 kHz f.sub.sym) K.sub.□ 24 (96 kHz BW) 80 (320 kHz BW) K.sub.□ 20 (K.sub.sub □ 480) 60 (K.sub.sub □ 4,800) W.sub.T 1.92 19.2 MHz f.sub.DAC 3.84 Msps IQ 38.4 Msps IQ N.sub.DAC 960 (2,880 bytes) 9,600 (28.125 KB) SF.sub.LPF 3 (960 kHz passband) 3 (9.6 MHz passband)

[0105] Assuming 12 bits per in-phase (I) and quadrature (Q) rail at the output of the quantizer 108, a single time symbol vector 113 requires transmission of 2,880 bytes (2.8125 KB) over the beacon communication bus 114 for the LMS Scenario, and 28.125 KB for the 2.4 GHz Ch. 13 Scenario. Assuming the time symbol vectors 113 are updated once per second, the NOC 403 requires the beacon communication bus 114 to support a 23.04 kbps link for the LMS scenario, and a 230.4 kbps link for the 2.4 GHz Ch. 13 Scenario. These rates are achievable in low-cost networks.

[0106] FIG. 10A-FIG. 10D show the frequency and time response of two 2.4 GHz Ch. 13 beacons with parameters consistent with Table 1. FIG. 10A shows the flat beacon frequency response over the 19.2 MHz signal passband, with negligible sidebands. FIG. 10B shows the close-in spectral distribution of the beacons—essentially, a line spectrum with 4 kHz line separation and flat line power—and illustrates the layout of the 320 kHz clusters within the signal. FIG. 10C and FIG. 10D show the instantaneous power of the beacon, and demonstrates the effect of optimized code phasings on the transmitted signal, one of the features of the disclosure. In FIG. 10C, the beacon is generated using inner subcarrier vector 104 and outer subcarrier vector 106 code phases with uniform-random phase distribution. This beacon has a PAPR of 11.2 dB, consistent with a bandlimited complex-Gaussian waveform. In FIG. 10D, the other beacon is generated using 80-element inner subcarrier vector 104 and 60-element outer subcarrier vector 106 code phases taken from 1,000-member code libraries 102, in which the phases are optimized to minimize kurtosis of their underlying time series. This beacon has a PAPR of 3.15 dB.

Exemplary Beacon Reception System

[0107] FIG. 11 shows a receiver system used in one aspect of the disclosure. When directed to collect a snapshot, e.g., based on prompts from the NOC 403 over a wireless communication transceiver 402, or at scheduled snapshot collection times, the receiver system performs reception operations 300 to generate a downconverted and sampled data stream covering the beacon transmission band. The reception operations shown in FIG. 11 can comprise optional bandpass filtering (BPF), followed by low-noise amplification (LNA), in-phase and quadrature (IQ) downconversion to complex baseband, represented as multiplication of the LNA output signal by a complex LO output signal, dual lowpass filtering (LPF), and dual analog-to-digital conversion (ADC) operations, resulting in a complex IQ sampled received data stream. However, the reception operations 300 can be implemented in many other ways, including two-stage superheterodyne reception operations 300 that downconvert the beacon transmission band to real-IF representation, and hybrid analog-digital reception operations 300 that downconvert a larger frequency band containing the beacon transmission band, and implements a digital drop receiver to generate a reduced-bandwidth, decimated signal centered on the beacon transmission band.

[0108] When needed at scheduled intervals, or given prompts from the NOC 403 over a wireless communication transceiver 402, the receiver system then performs a data snapshot collection 301 operation, which generates a snapshot 302 comprising the data provided by the reception operations 300 at a reception time and over a snapshot 302 time duration, shown in FIG. 11 as time center t.sub.R.sup.ref, and time duration T.sub.R also included as a time-stamp with that snapshot 302. In aspects with slotted beacon formats, a prefix and suffix with duration T.sub.prefix and T.sub.suffix, respectively, are also collected as part of the snapshot 302, in order to encompass inter-slot interference introduced by timing offset between the beacon transmitter and receiver. The snapshot 302 is then sent to a position/timing (P/T) solution generator 303 over a snapshot communication bus 304. In the network shown in FIG. 1, the snapshot communication bus 304 can be connected to the NOC 403 over a wireless communication transceiver 402, or it can connect to a P/T solution generator 303 on-board the user 400. In the networks shown in FIG. 2 and FIG. 3, the snapshot communication bus 304 is connected to the NOC 403 over a link similar to links supporting the beacon communication bus 114. However, the snapshot communication bus 304 typically requires a higher-rate link than the beacon communication bus 114.

[0109] In general, the receiver LO(s) and ADC samplers used in the reception operation 300 are locked to a system clock 214 with rate offset ε.sub.R and timing offset τ.sub.R.sup.ref=t.sub.R.sup.ref−t.sub.UTC(t.sub.R.sup.ref), unique for each receiver, where t.sub.R.sup.ref is the receiver time estimate at UTC time t.sub.UTC(t.sub.R.sup.ref).

[0110] In some aspects of the disclosure, the NOC 403 also provides synchronization data 216 that can be used to bring the receiver system clock 214 into synchronization for subsequent time-stamped snapshots 302. In other aspects, the receiver obtains coarse synchronization information from the NOC 403 over a wireless communication transceiver 402. In additional aspects of the disclosure, the receiver performs coarse synchronization operations to determine the approximate center frequency and (for slotted beacon formats) slot transition time of the beacons, prior to the snapshot collection 301. The coarse frequency and timing information can then be used to adjust the timing carrier offset of the ADC output signal, or the receiver clock driving the LO and ADC; or simply conveyed to the NOC 403, along with time-stamped snapshot 302. In this last case, the frequency and timing offset is included in the time-stamped snapshot 302, for use by the P/T solution generator 303. The receiver may stream data to the P/T solution generator 303 (which may be remote or on the receiver platform itself), or may sparsely capture time-stamped snapshot 302 of frequency-and-timing aligned data, e.g., at the start of processing or as required/requested by the NOC 403.

Cold-Start Reception and Resilient Geo-Observable Estimation Operations

[0111] FIG. 12 shows exemplary cold-start operations used at the P/T solution generator 303 to resiliently detect and estimate geo-observables of co-channel beacons contained within individual received snapshots 302, applicable to the network shown in FIG. 1. The snapshot 302 is first processed to estimate the FOA centroid 612, described in FIG. 15, i.e., the center of the beacon FOA's observed in the downconverted snapshot 303, which is roughly offset by the error between the target and actual LO in the reception operations 300. This can be accomplished in a number of manners; for example, using spectral analysis tools to determine the center of the received beacon spectrum, or to exploit the nearly-rectangular shape of the beacon signals, e.g., by detecting and estimating the up-edge and down-edge of the signal spectrum. The chief goal of this step is to determine the FOA centroid 612 to within a fraction of a subcarrier, e.g., ±500 Hz (±1/8 subcarrier) for the beacon transmission parameters given in Table 1. In some aspects with slotted beacon formats, the snapshot 303 is optionally analyzed to estimate a timing offset estimate {circumflex over (n)}.sub.R, e.g., by performing cross-slot time correlation operations to detect the slot transition.

[0112] The snapshot 303 is then frequency-shifted to remove the estimated FOA centroid 612, and if needed time-shifted to remove the estimated timing offset, and channelized into subcarriers and time symbols covering the active snapshot bandwidth and duration 502, creating a channelized snapshot 553, described in FIG. 13, defined over subcarriers and time symbols. The channelized snapshot 553 is then stacked into an N.sub.data×M.sub.DoF windowed data matrix 601, described in FIG. 15, where M.sub.DoF is the degrees-of-freedom (DoF) of the data matrix and N.sub.data is the number of data samples in the matrix, and where the beacons are redundant in the DoF dimension; and whitened over the DoF dimension, e.g., using a QR decomposition (QRD) operation 503, thereby creating a whitened snapshot matrix 604, described in FIG. 15. In some aspects, FOA clutter and FOA centroid computation operations 504 are performed to generate a FOA clutter spectrum 610, described in FIG. 15, and used to improve detection of the received beacons, and optionally to improve the estimate of the FOA centroid 612. In this case, the improved FOA centroid 612 estimate can also be used to regenerate 505 the channelized snapshot 553 and the whitened snapshot matrix 604.

[0113] The whitened snapshot matrix 604 is then passed through an FFT/IFFT mechanized resilient least-square (LS) search operation 507 to form a least-squares (LS) TOA-FOA surface 613, described in FIG. 15, which provides a detection statistic as a function of candidate TOA and FOA values, and to search the LS TOA-FOA surface 613 to detect each beacon in the whitened snapshot matrix 604, using known data dimension codes 506 applied during beacon time-symbol generation operations shown in FIG. 7, and to determine the observed TOA and FOA, and LS SINR (estimation quality) of the detected beacons. The maximizing TOA and FOA, and the maximal LS SINR are then refined 508, for example, using polynomial fit to the surface peak or parametric search operations in the vicinity of the detected TOA-FOA surface peak. Optionally, a copy-aided parametric estimation method 510 is used to further refine TOA and FOA geo-observables and resolve detector ambiguities 510, using known DoF dimension codes 506 applied during beacon generation operations. In some aspects, results of the copy-aided parametric estimation method 510 is used to further improve estimates of the FOA centroid 612, and optionally the estimate of the timing offset, and the channelized snapshot 553 is regenerated 505 using the improved FOA centroid 612 and optional timing offset estimates, to further improve the geo-observable estimates.

[0114] The DoF and data dimensions, and the data dimension codes 506 and DoF dimension codes 509 are set based on the particular form of redundancy exploited in the FFT/IFFT mechanized resilient least-square (LS) search operation 507, as determined by the stacking operation performed in the data stacking and whitening operation 503. For example, in one aspect where the channelized snapshot is stacked over the inner-code dimension, data dimension N.sub.data=K.sub.1N.sub.sym, and the data dimension codes 506 are the K.sub.1×1 outer subcarrier vector 106 phases ; and DoF dimension M.sub.DoF=K.sub.0 and the DoF dimension codes 509 are the K.sub.0×1 inner subcarrier vector 104 phases . In a second aspect where the channelized snapshot is stacked over the outer-code dimension, data dimension N.sub.data=K.sub.0N.sub.sym, and the data dimension codes 506 are the K.sub.0×1 inner subcarrier vector 104 phase ; and DoF dimension M.sub.DoF=K.sub.1 and the DoF dimension codes 509 are the K.sub.1×1 outer subcarrier vector 106 phases .

[0115] The key components of this procedure are described in more detail in the next subsections.

Channelization Operations

[0116] FIG. 13 depicts exemplary operations used to implement the FOA centroid 612 and optional timing estimate removal and channelization operations 502 in one aspect. The system first performs a FOA centroid and optional timing offset removal operation 550, to remove the estimated FOA centroid 612 {circumflex over (α)}.sub.R and optional estimated timing offset {circumflex over (n)}.sub.R from the snapshot 303. In some aspects with slotted beacon formats, the prefix and suffix are also discarded prior during the FOA centroid and optional timing offset removal operation 550. The snapshot 303 then has duration T.sub.R=N.sub.repT.sub.sym=N.sub.RT.sub.ADC, where N.sub.rep and N.sub.R are integers and T.sub.ADC=1/f.sub.ADC is the ADC sampling period in the receiver's frame of reference.

[0117] The snapshot 303 is then separated into frequency subcarriers and time symbols 551 using a sparse, overlapped, optionally frequency-offset, windowed DFT with overlap time T.sub.sym N.sub.ADCT.sub.ADC, sparsity factor Q.sub.sym, DFT length N.sub.DFT=Q.sub.symN.sub.ADC, and channelizer window {w.sub.R(m.sub.ADC)}.sub.m.sub.ADC.sub.=0.sup.N.sup.DFT.sup.−1, and with frequency offset factor f.sub.sym/2 if K.sub.sub is even. The DFT output bins corresponding to the active beacon subcarrier frequencies are then selected, and a snapshot equalizer operation 552 is applied to those bins. These operations generates a channelized snapshot 553 x.sub.sub (k.sub.sub,n.sub.sym) given by

[00008]$\begin{array}{cc}{x}_{sub}\left(k_{sub},n_{sym}\right)=G_R\left(k_{sub},n_{sym};{\hat{\alpha }}_R\right){\tilde{x}}_{sub}\left(k_{sub},n_{sym}\right),\{\begin{array}{c}{k}_{sub}=0,.Math.,K_{sub}-1\\ n_{sym}=0,.Math.,N_{sym}-1,\end{array}& \left(\mathrm{Eq6}\right)\end{array}$$\begin{array}{cc}{\tilde{x}}_{sub}\left(k_{sub},n_{sym}\right)=\stackrel{N_{\mathrm{DFT}}-1}{\underset{m_{\mathrm{ADC}}=0}{.Math.}}{w}_R\left(m_{ADC}\right)x_R\left(t_R\left(n_{sym}\right)+m_{ADC}{T}_{ADC}\right)\times \exp \left\{-j2\pi \left(f\left(k_{sub}\right)+{\hat{\alpha }}_R\right)\left(t_R\left(n_{sym}\right)+m_{ADC}{T}_{ADC}\right)\right\},& \left(\mathrm{Eq7}\right)\end{array}$

where N.sub.sym=N.sub.rep−Q.sub.sym+1 is the number of time symbols in the channelized snapshot and

[00009]$t_R\left(n_{sym}\right)=t_R^{\mathrm{ref}}+\left(n-\frac{N_{sym}-1}{2}\right)T_{sym}$

is the symbol n.sub.sym DFT time-center in the receiver's field of reference, and where G.sub.R(k.sub.sub,n.sub.sym;{circumflex over (α)}.sub.R) are snapshot channelizer equalizer 552 weights that remove effects caused by at least the carrier operation 552, and optionally filtering effects of the reception operations 300. Preferentially, the snapshot equalizer 552 weights are given by

[00010]$\begin{array}{cc}{G}_R\left(k_{sub},n_{sym};{\hat{\alpha }}_R\right)=H_R^{-1}\left(\left(1-\frac{{\hat{\alpha }}_R}{f_T}\right)\left(f\left(k_{sub}\right)+{\hat{\alpha }}_R\right)\right){\delta }_R\left(k_{sub},n_{sym},{\hat{\alpha }}_R\right),& \left(\mathrm{Eq8}\right)\end{array}$

where H.sub.R(f) is the aggregate frequency response of the reception operations 300, and where δ.sub.R(k.sub.sub,n.sub.sym;{circumflex over (α)}.sub.R) removes dispersive effects of the FOA centroid 612 removal operation 502,

[00011]$\begin{array}{cc}{\delta }_R\left(k_{sub},n_{sym};{\hat{\alpha }}_R\right)=\exp \left\{-j2\pi f\left(k_{sub}\right)\left(t_R\left(n_{sym}\right)-t_R^{\mathrm{ref}}\right)\frac{{\hat{\alpha }}_R}{f_T}\right\}.& \left(\mathrm{Eq9}\right)\end{array}$

The aggregate frequency response term is optional, and can be based on modeling of the reception operations 300; or derived from calibration operations performed by the receiver or network, for example, as described in Pattabiraman 2012, and can be based on the magnitude or complex value of those reception operations 300.

[0118] Table 2 lists receiver and channelizer parameters compatible with the beacon generation and transmission parameters shown in Table 1. The receiver assumes a dual-ADC sampling rate of 3.84 million samples per second (Msps) for the LMS Scenario, and 30.72 Msps for the 2.4 GHz Ch. 13 Scenario, with sufficient antialiasing filtering to provide a 2 MHz and 20 MHz protected two-way passband, respectively, covering the active bandwidth of the beacons with ±40 kHz and ±400 kHz of guard band for LO uncertainty, respectively. A mixed-radix DFT with factor-of-four sparsity (Q.sub.sym=4) is assumed in both scenarios, and a separation of 250 μs between successive DFT's. The Table further assumes a 10 millisecond snapshot encompassing 16 symbol repetitions, 13 of which are used in subsequent geo-observable estimation operations, for both scenarios.

TABLE-US-00002 TABLE 2 Exemplary Receiver, Channelizer Parameters Parameter LMS Scenario 2.4 GHz Ch. 13 Scenario SF.sub.LPF 2.84 (1 MHz passband) 2.072 (10 MHz passband) f.sub.ADC 3.84 Msps 30.72 Msps W.sub.R 2 MHz 20 MHz W.sub.guard 40 kHz 400 kHz T.sub.sym 250 μs (N.sub.ADC □ 960) 250 μs (N.sub.ADC □ 9,600) Q.sub.sym 4 (480 subcarriers) 4 (4,800 subcarriers) N.sub.sym 13 (16 repetitions) 13 (repetitions) T.sub.R 10 ms (N.sub.R = 38,400) 10 ms (N.sub.R = 307,200) Snapshot size 112.5 KB 900 KB

[0119] Assuming dual-ADC precision of 12 bits per I and Q rail, consistent with a low-cost receiver front-end, the size of each snapshot is 112.5 KB for the LMS Scenario, and 900 KB for the 2.4 GHz Ch. 13 Scenario. Assuming a snapshot is collected 303 once per second, backhaul 303 of ADC output data to the P/T solution generator requires a snapshot communication bus 304 that can support a 0.922 Mbps one-way data-rate for the LMS Scenario, and a 7.37 Mbps one-way data-rate for the 2.4 GHz Ch. 13 Scenario, well within capabilities of 4G cellular or 802.11 WLAN standards if the P/T solution generator 303 is in the NOC 403. Continuous backhaul of snapshots 302 to the P/T solution generator 303 over the snapshot communication bus 304 would require a factor of 100 higher data-rate, e.g., 92.2 Mbps and 737 Mbps, respectively, for the two scenarios, easily accomplished over Gbps Ethernet if the P/T solution generator 303 is on-board the user 400.

[0120] The FOA centroid 612 and optional timing estimate removal and channelization operations 502 can be performed in a number of different manners, for example, using polyphase filtering methods, discrete filter banks centered on each subcarrier frequency, mixtures of radix-2 and non-radix-2 fast Fourier transform (FFT) and inverse-FFT methods, and so on.

[0121] FIG. 14 shows the frequency response of an exemplary DFT window developed for the LMS Scenario using an interpolated Parks-McClellan algorithm. The window passband is set to 1 kHz, corresponding ±136 Hz FOA offset expected in this band at the maximum class-1 sUAS airspeed of 44.7 meters/second (100 miles/hour), with ±864 Hz carrier offset margin after FOA centroid 612 removal operations. The window stopband is set to 3 kHz, sufficient to reject any inter-subcarrier interference due to FOA offset from the receiver center frequency. As shown in this FIG., the window provides 3.7 dB rejection within its design passband (much less within ±136 Hz), and over 58 dB of rejection within its design stopband. The DFT window designed for the 2.4 GHz Ch. 13 Scenario has nearly identical performance.

[0122] Assuming synchronized beacon transmitters in the network shown in FIG. 1, such that (t.sub.UTC)≡t.sub.UTC and defining t.sub.UTC.sup.reft.sub.UTC(t.sub.R.sup.ref) as the actual UTC time (to be estimated in positioning/timing algorithms) at receiver reference time t.sub.R.sup.ref, and further assuming short snapshots 302 are collected 301, then the channel link gain between the beacon transmitter deployed at network node 401 and a receiver deployed at a user 400 is approximated by (t.sub.UTC)≈(r.sub.UTC.sup.ref), and the TOA and FOA of the beacon received at that user 400 is approximated by

[00012]$\begin{array}{cc}{\tau }_{T\left(\ell \right)R}\left(t_{UTC}\right)\approx {\tau }_{T\left(\ell \right)R}\left(t_{UTC}^{\mathrm{ref}}\right)+\left(t_{UTC}-t_{UTC}^{\mathrm{ref}}\right){\epsilon }_{T\left(\ell \right)R}\left(t_{UTC}^{\mathrm{ref}}\right),& \left(\mathrm{Eq10}\right)\end{array}$$\begin{array}{cc}{\alpha }_{T\left(\ell \right)R}\left(t_{UTC}\right)\approx -f_T{\epsilon }_{T\left(\ell \right)R}\left(t_{UTC}^{\mathrm{ref}}\right),\mathrm{where}{\tau }_{T\left(\ell \right)R}\left(t_{UTC}^{ref}\right)=\frac{1}{c}{.Math.p_{T\left(\ell \right)R}\left(t_{UTC}^{\mathrm{ref}}\right).Math.}_2^2\mathrm{and}{\epsilon }_{T\left(\ell \right)R}\left(t_{UTC}^{\mathrm{ref}}\right)=-\frac{1}{c}{u}_{T\left(\ell \right)R}^T\left(t_{UTC}^{\mathrm{ref}}\right)v_R\left(t_{UTC}^{\mathrm{ref}}\right)& \left(\mathrm{Eq11}\right)\end{array}$

are the TOA and differential TOA (DTOA) of beacon at time t.sub.UTC.sup.ref, and where (t.sub.UTC)=−p.sub.R(t.sub.UTC) is the observed position of network node 401 at the user 400, and (t.sub.UTC)=(t.sub.UTC)/∥(t.sub.UTC)∥.sub.2 is the observed line-of-bearing (LOB) from the user 400 to network node 401 . The TOA, DTOA, and FOA observed at the ADC sampler in the user 400 reception operations 300 are then given by

[00013]$\begin{array}{cc}{\tilde{\tau }}_{T\left(\ell \right)R}\left(t_R\right)={\tilde{\tau }}_{T\left(\ell \right)R}\left(t_R^{\mathrm{ref}}\right)+\left(t_R-t_R^{\mathrm{ref}}\right){\tilde{\epsilon }}_{T\left(\ell \right)R}\left(t_R^{\mathrm{ref}}\right),& \left(\mathrm{Eq12}\right)\end{array}$$\begin{array}{cc}{\tilde{\tau }}_{T\left(\ell \right)R}\left(t_R^{ref}\right)\approx {\tau }_{T\left(\ell \right)R}\left(t_{UTC}^{\mathrm{ref}}\right)+\left(t_R^{\mathrm{ref}}-t_{UTC}^{\mathrm{ref}}\right),& \left(\mathrm{Eq1}3\right)\end{array}$$\begin{array}{cc}{\tilde{\epsilon }}_{T\left(\ell \right)R}\left(t_R^{\mathrm{ref}}\right)=\frac{{\epsilon }_R}{1+{\epsilon }_R}+\frac{1}{1+{\epsilon }_R}{\epsilon }_{T\left(\ell \right)R}\left(t_{UTC}^{\mathrm{ref}}\right),& \left(\mathrm{Eq1}4\right)\end{array}$$\begin{array}{cc}{\tilde{\alpha }}_{T\left(\ell \right)R}\left(t_R^{\mathrm{ref}}\right)=-f_T{\tilde{\epsilon }}_{T\left(\ell \right)R}\left(t_R^{\mathrm{ref}}\right).& \left(\mathrm{Eq1}5\right)\end{array}$

[0123] Further assuming accurate equalization of user 400 reception operations 300, and ideal suppression of inter-subcarrier interference by the channelizer window used in the FOA centroid 612 and optional timing estimate removal and channelization operations 502, channelized snapshot 553 x.sub.sub(k.sub.sub,n.sub.sym) is approximated by

[00014]$\begin{array}{cc}{x}_{sub}\left(k_{sub},n_{sym}\right)\approx i_{sub}\left(k_{sub},n_{sym}\right)+\underset{\ell =1}{\overset{L}{.Math.}}{s}_{sub}\left(k_{sub},n_{sym};\ell \right),& \left(\mathrm{Eq16}\right)\end{array}$$\begin{array}{cc}{s}_{sub}\left(k_{sub},n_{sym};\ell \right)\approx a_{sub}\left(\ell \right)d_{sub}\left(k_{sub},n_{sym};{\tilde{\tau }}_{T\left(\ell \right)R}\left(t_R^{\mathrm{ref}}\right),{\tilde{\alpha }}_{T\left(\ell \right)R}\left(t_R^{\mathrm{ref}}\right)-{\hat{\alpha }}_R;\ell \right),& \left(\mathrm{Eq1}7\right)\end{array}$

where α.sub.sub() is the end-to-end beacon channelizer output gain,

[00015]$\begin{array}{cc}{a}_{sub}\left(\ell \right)\approx W_R\left(e^{-j2\pi \left({\tilde{\alpha }}_{T\left(\ell \right)R}\left(t_R^{\mathrm{ref}}\right)-{\hat{\alpha }}_R\right)T_{ADC}}\right){\tilde{a}}_{T\left(\ell \right)R}){\tilde{a}}_{T\left(\ell \right)R}\left(t_R^{\mathrm{ref}}\right),& \left(\mathrm{Eq18}\right)\end{array}$$\begin{array}{cc}{\tilde{a}}_{T\left(\ell \right)R}\left(t_R^{\mathrm{ref}}\right)=g_{T\left(\ell \right)R}\left(t_{UTC}^{\mathrm{ref}}\right)\sqrt{\frac{P_{T\left(f\right)}}{K_{sub}}}\exp \left\{j\left({\phi }_{T\left(\ell \right)}-{\phi }_R-2\pi f_T{\tilde{\tau }}_{T\left(\ell \right)R}\left(t_R^{\mathrm{ref}}\right)\right)\right\},& \left(\mathrm{Eq1}9\right)\end{array}$

and where d.sub.sub(k.sub.sub,n.sub.sym;τ,α;) is the network node 401 beacon (“beacon ”) at candidate observed TOA τ and observed FOA α,

[00016]$\begin{array}{cc}{d}_{sub}\left(k_{sub},n_{sym};\tau ,a;\ell \right)=\exp \left\{j\left({\theta }_{T\left(\ell \right)}\left(k_{sub}\right)-2\pi \left(f\left(k_{sub}\right)\tau -\left(t_R\left(n_{\mathrm{sym}}\right)-t_R^{\mathrm{ref}}\right)a\right)\right)\right\}\times \exp \left\{j2\pi f\left(k_{sub}\right)\left(t_R\left(n_{sym}\right)-t_R^{\mathrm{ref}}\right)\frac{a}{f_T}\right\},\mathrm{and}{W}_R\left(e^{j2\pi f}\right)=\underset{m}{.Math.}{w}_R\left(m\right)\exp \left\{-j2\pi \mathrm{mf}\right\}& \left(\mathrm{Eq20}\right)\end{array}$

is the analytic discrete Fourier transform of the channelizer window. Assuming additive white Gaussian noise (AWGN) with noise density N.sub.0 at the LPF input and ideal LPF equalization, the background interference i.sub.sub(k.sub.sub,n.sub.sym) has identical power R.sub.i.sub.sub.sub.i.sub.sub≈f.sub.ADCN.sub.0∥w.sub.R(m.sub.ADC)∥.sub.2.sup.2 on each channelizer output subcarrier, and beacon has identical SNR ≈|a.sub.sub()|.sup.2/R.sub.i.sub.sub.sub.i.sub.sub on each channelizer output subcarrier.

[0124] Using the channelized snapshot 553 model given in (Eq16)-(Eq20), the observed geo-observables {{tilde over (τ)}.sub.T(.sub.)R(t.sub.R.sup.ref),(t.sub.R.sup.ref)} can in principle be estimated by correlating channelized snapshot 553 x.sub.sub(k.sub.sub,n.sub.sym) against {d.sub.sub(k.sub.sub,n.sub.sym;). Moreover, if

[00017]$f_T\frac{\pi }{2}{N}_{sym}{K}_{sub}.Math.\alpha .Math.,$

then the second dispersive term in (Eq20) can be ignored, and this correlation can be efficiently mechanized using DFT and inverse-DFT (IDFT) methods. However, at high receive SNR this correlation will yield a poor result, or will require a high time-bandwidth product N.sub.symK.sub.sub to remove cross-correlation between co-channel beacons. This problem can be overcome by exploiting the spectral or temporal redundancy imposed in this signal at the transmitter, to implement resilient TOA-FOA estimators. This procedure is described below.

Exemplary Resilient Fine Least-Squares TOA-FOA Estimator

[0125] FIG. 15 depicts one aspect employing resilient TOA-FOA estimation, referred to here as the fine least-squares procedure. The channelizer snapshot 553 is first stacked along it's inner code dimension 600, generating K.sub.0×1 vector signal x.sub.0(k.sub.1,n.sub.sym)=[x.sub.sub(K.sub.0k.sub.1+k.sub.0,n.sub.sym)].sub.k.sub.0.sub.=0.sup.K.sup.0.sup.−1, enumerated over K.sub.1 outer subcarrier channel indices k.sub.1=0, K.sub.1−1 and N.sub.sym time symbol indices n.sub.sym=0, N.sub.sym−1. Using (Eq16)-(Eq20), this signal can be modeled as

[00018]$\begin{array}{cc}{x}_0\left(k_1,n_{sym}\right)\approx i_0\left(k_1,n_{\mathrm{sym}}\right)+\underset{\ell =1}{\overset{L}{.Math.}}{s}_0\left(k_1,n_{\mathrm{sym}};\ell \right),\mathrm{where}{\left\{s_0\left(k_1,n_{\mathrm{sym}};\ell \right)\right\}}_{\ell =1}^L& \left(\mathrm{Eq2}1\right)\end{array}$

are the K.sub.0×1 inner-code stacked beacon signal vectors, respectively,

[00019]$\begin{array}{cccc}{s}_0\left(k_1,n_{sym};\ell \right)& =& \left({\delta }_0\left(n_{sym};{\tilde{\alpha }}_{T\left(\ell \right)R}-{\hat{\alpha }}_R\right)\circ a_0\left({\tilde{\tau }}_{T\left(\ell \right)R};\ell \right)\right)& \\ & \times & \left({\delta }_1\left(k_1,n_{sym};{\tilde{\alpha }}_{T\left(\ell \right)R}-{\hat{\alpha }}_R\right)d_1\left(k_1,n_{sym};{\tilde{\tau }}_{T\left(\ell \right)R},{\tilde{\alpha }}_{T\left(\ell \right)R}-{\hat{\alpha }}_R;\ell \right)\right),& \left(\mathrm{Eq22}\right)\\ & \approx & a_0\left({\tilde{\tau }}_{T\left(\ell \right)R};\ell \right)d_1\left(k_1,n_{sym};{\tilde{\tau }}_{T\left(\ell \right)R},{\tilde{\alpha }}_{T\left(\ell \right)R}-{\hat{\alpha }}_R;\ell \right),& \left(\mathrm{Eq23}\right)\end{array}$

and where a.sub.0(τ;) and d.sub.1(k.sub.1,n.sub.sym;τ,α;) is the K.sub.0×1 beacon outer-code spectral signature at trial TOA τ and scalar inner-code signal at trial TOA τ and FOA α, respectively,

a.sub.0(τ;)=a.sub.sub()[exp{j((k.sub.0)−2πf.sub.0(k.sub.0)τ)}].sub.k.sub.1.sub.=0.sup.K.sup.1.sup.−1,  (Eq24)

d.sub.1(k.sub.1,n.sub.sym;τ,α;=exp{j((k.sub.0)−2π(f.sub.1(k.sub.1)τ−(t.sub.R(n.sub.sym)−t.sub.R.sup.ref)α))},  (Eq25)

and the dispersive terms δ.sub.0(n.sub.sym;α) and δ.sub.1(k.sub.1,n.sub.sym;α) are given by

[00020]$\begin{array}{cc}{\delta }_0\left(n_{sym};\alpha \right)={\left[\exp \left\{j2\pi f_0\left(k_0\right)\left(t_R\left(n_{sym}\right)-t_R^{\mathrm{ref}}\right)\frac{\alpha }{f_T}\right\}\right]}_{k_0=0}^{K_0-1},& \left(\mathrm{Eq26}\right)\end{array}$$\begin{array}{cc}{\delta }_1\left(k_1,n_{sym};\alpha \right)=\exp \left\{j2\pi f_1\left(k_1\right)\left(t_R\left(n_{sym}\right)-t_R^{\mathrm{ref}}\right)\frac{a}{f_T}\right\},& \left(\mathrm{Eq27}\right)\end{array}$

respectively. If

[00021]$f_T\frac{\pi }{2}{K}_0{N}_{sym}.Math.{\tilde{\alpha }}_{T\left(\ell \right)R}-{\hat{\alpha }}_R.Math.,$

then δ.sub.0(n.sub.sym;−{circumflex over (α)}.sub.R)≈1.sub.K.sub.0 and s.sub.0(k.sub.1,n.sub.sym;) is nondispersive over the inner-code dimension. Similarly, if

[00022]$f_T\frac{\pi }{2}{K}_{sub}{N}_{sym}.Math.{\tilde{\alpha }}_{T\left(\ell \right)R}-{\hat{\alpha }}_R.Math.,$

then δ.sub.1(k.sub.1,n.sub.sym;α)≈1 and (Eq23) holds closely. Assuming AWGN background noise and ideal LPF equalization, the K.sub.0×1 background interference vector i.sub.0(k.sub.1,n.sub.sym)=[i.sub.sub(K.sub.0k.sub.1+k.sub.0,n.sub.sym)].sub.k.sub.0.sub.=0.sup.K.sup.0.sup.−1 has asymptotic autocorrelation matrix (ACM) R.sub.i.sub.0.sub.i.sub.0.fwdarw.R.sub.i.sub.sub.sub.i.sub.subI.sub.K.sub.0 on each inner subcarrier channel.

[0126] This model is closely analogous to a multi-element antenna array with M.sub.DoF degrees-of-freedom (DoF's), where M.sub.DoF=K.sub.0 is the stacking or “DoF” dimension. Similar to an array, it allows the outer-code signals to be detected and separated with an output (despread) SINR approximated by ˜(M.sub.DoF−L+1 for in presence of strong interference from co-channel beacons ( 1, where ′≠),using well-known, mature linear signal separation methods, e.g., least-squares (LS) algorithms, referred to as code nulling in Agee 2000. The method sacrifices one despreader DoF to null each strong signal in the environment, and uses the despreader's remaining DoF's to improve the output SNR of the intended signal. It also admits superresolution geo-observable estimators with accuracy that scales with this output SINR.

[0127] The inner-code stacked signal is then formed into K.sub.1N.sub.sym×K.sub.0 windowed data matrix 601 X.sub.0=[√{square root over (w.sub.TOA(k.sub.1)w.sub.FOA(n.sub.sym))}x.sub.0.sup.T(k.sub.1,n.sub.sym)], where tapering windows 602 {w.sub.TOA(k.sub.1)}.sub.k.sub.1.sub.=0.sup.K.sup.1.sup.−1 and {w.sub.FOA(n.sub.sym)}.sub.n.sub.sym.sub.=0.sup.N.sup.sym.sup.−1 satisfy Σw.sub.TOA(k.sub.1)=Σw.sub.FOA(n.sub.sym)=1, and passed through a whitening operation 603, to generate whitened data matrix 604 Q.sub.0=[q.sub.0.sup.T(k.sub.1,n.sub.sym)]. The whitening operation 603 can be performed, for example, using QR decomposition (QRD) {Q.sub.0,R.sub.0}=QRD(X.sub.0), given by Q.sub.0=X.sub.0R.sub.0.sup.−1, where R.sub.0 is an upper-triangular matrix given by R.sub.0=chol{X.sub.0.sup.HX.sub.0}, and where chol{∘} is the Cholesky factor operation. The QRD can be implemented using a number of efficient methods known to those of ordinary skill in the arts, e.g., modified Gram-Schmidt orthogonalization (MGSO). Other whitening operations, e.g., singular value decomposition, can also be used to whiten X.sub.0. The whitened data matrix 604 is then used to form SINR-revealing FLS TOA-FOA surface 613 {circumflex over (γ)}.sub.FLS(τ,α;), by computing intermediate K.sub.0×1 FLS FOA vector 606

[00023]$\begin{array}{cc}{r}_{FLS}\left(k_1;a\right)=\sqrt{w_{TOA}\left(k_1\right)}\underset{n_{sym}=0}{\overset{N_{sym}-1}{.Math.}}\sqrt{w_{FOA}\left(n_{sym}\right)}{q}_1^{*}\left(k_1,n_{sym}\right)\exp \left\{j2\pi n_{sym}{T}_{sym}a\right\},& \left(\mathrm{Eq28}\right)\end{array}$

on each inner subcarrier channel using an DFT bank 605; computing K.sub.0×1 whitened FLS linear combiner vector

[00024]$\begin{array}{cc}{\hat{u}}_{\mathrm{FLS}}\left(\tau ,\alpha ;\ell \right)=\underset{k_1=0}{\overset{K_1-1}{.Math.}}{r}_{FLS}\left(k_1;\alpha \right)\exp \left\{j\left({\theta }_{T_1\left(\ell \right)}\left(k_1\right)-2\pi k_1{f}_{\mathrm{sym}}{K}_0\tau \right)\right\},& \left(\mathrm{Eq29}\right)\end{array}$

for each candidate FOA and data-dimension code 506, for this surface the outer-code, using an IDFT bank; and computing SINR-revealing FLS TOA-FOA surface 613

[00025]$\begin{array}{cc}{\hat{\gamma }}_{\mathrm{FLS}}\left(\tau ,\alpha ;\ell \right)=\frac{{.Math.{\hat{u}}_{FLS}\left(\tau ,\alpha ;\ell \right).Math.}_2^2}{1-{.Math.{\hat{u}}_{FLS}\left(\tau ,\alpha ;\ell \right).Math.}_2^2}.& \left(\mathrm{Eq30}\right)\end{array}$

Optionally, the FLS FOA vector 606 also yields FLS FOA clutter spectrum 610 {circumflex over (γ)}.sub.FLS-clutter(α)={circumflex over (η)}.sub.FLS-clutter (α)/(1−{circumflex over (η)}.sub.FLS-clutter(α)), where

[00026]$\begin{array}{cc}{\hat{\eta }}_{FLS-cl\mathrm{utter}}\left(\alpha \right)=\stackrel{K_1-1}{\underset{k_1=0}{.Math.}}{.Math.{\hat{r}}_{FLS}\left(k_{sub},\alpha \right).Math.}_2^2,& \left(\mathrm{Eq31}\right)\end{array}$

which can be used to compute FLS deflection statistic {circumflex over (d)}.sub.FLS(τ,α;)={circumflex over (γ)}.sub.FLS(τ,α;)/{circumflex over (γ)}.sub.FLS-clutter(α), a particularly useful statistic in TOA-FOA spectra containing multiple significant peaks, e.g., due to specular multipath. In some aspects, the clutter statistic is also used to improve the FOA centroid 612, e.g., using formula

[00027]$\begin{array}{cc}{\hat{\alpha }}_R\leftarrow {\hat{\alpha }}_R+\frac{\underset{\alpha }{.Math.}(x{\hat{\gamma }}_{FLS-\mathrm{clutter}}\left(\alpha \right)}{\underset{\alpha }{.Math.}{\hat{\gamma }}_{FLS-\mathrm{clutter}}\left(\alpha \right)},& \left(\mathrm{Eq32}\right)\end{array}$

which can be used to regenerate the channelized snapshot 553.

[0128] In absence of substantive multipath, the TOA-FOA estimate is then given by

[00028]$\begin{array}{cc}\left({\hat{\tau }}_{T\left(\ell \right)R},{\hat{\alpha }}_{T\left(\ell \right)R}-{\hat{\alpha }}_R\right)=\arg \underset{\begin{array}{c}0\le \tau <T_{sym}/K_0\\ .Math.\alpha .Math.\le f_{\mathrm{sym}}\end{array}}{\max }{\hat{\gamma }}_{FLS}\left(\tau ,\alpha ;\ell \right),& \left(\mathrm{Eq33}\right)\end{array}$

and the maximal TOA-FOA surface 613 value is a metric of the SINR of the FLS combiner output signal at the estimated TOA and FOA, {circumflex over (γ)}.sub.FLS()={circumflex over (γ)}.sub.FLS(−{circumflex over (α)}.sub.R; ). The inner-stacked spectral signature is optionally estimated by â.sub.0()=(C.sub.0.sup.Hû.sub.FLS(,−{circumflex over (α)}.sub.R;)) where C.sub.0=R.sub.0.sup.−1. The TOA and FOA error variances are further optionally estimated by {circumflex over (σ)}.sub.(∘).sup.2()/{circumflex over (γ)}.sub.FLS(), where

[00029]$\begin{array}{cc}{\sigma }_{\mathrm{TOA}}^2=\frac{{.Math.w_{\mathrm{TOA}}\left(k_1\right).Math.}_2^2{.Math.w_{\mathrm{FOA}}\left(n_{\mathrm{sym}}\right).Math.}_2^2}{8{\pi }^2\left(.Math.f_1^2\left(k_1\right).Math.-{.Math.f_1\left(k_1\right).Math.}^2\right)}\ge \frac{3}{4{\pi }^2}\frac{T_{\mathrm{sym}}^2}{K_1{N}_{\mathrm{sym}}{K}_0^2\left(K_1^2-1\right)},& \left(\mathrm{Eq34}\right)\end{array}$$\begin{array}{cc}{\sigma }_{\mathrm{FOA}}^2=\frac{{.Math.w_{\mathrm{TOA}}\left(k_1\right).Math.}_2^2{.Math.w_{\mathrm{FOA}}\left(n_{\mathrm{sym}}\right).Math.}_2^2}{8{\pi }^2\left(.Math.t_R^2\left(n_{\mathrm{sym}}\right).Math.-{.Math.t_R\left(n_{\mathrm{sym}}\right).Math.}^2\right)}\ge \frac{3}{4{\pi }^2}\frac{f_{\mathrm{sym}}^2}{K_1{N}_{\mathrm{sym}}\left(N_{\mathrm{sym}}^2-1\right)},& \left(\mathrm{Eq35}\right)\end{array}$

and where the lower bounds in (Eq34)-(Eq35) are achieved for flat tapering windows 602. Moreover, the background TOA-FOA surface 613 values are a factor of 2∥w.sub.TOA(k.sub.1)∥.sub.2.sup.2∥w.sub.FOA((n.sub.sym)∥.sub.2.sup.2≥2/K.sub.1N.sub.sym below {circumflex over (γ)}.sub.FLS(), where the lower bound is also achieved for flat tapering windows 602. For this reason, flat tapering windows 602 are recommended in absence of channel multipath, and shaped tapering windows 602 recommended if multiple TOA-FOA surface 613 peaks are expected, e.g., due to strong specular multipath.

[0129] The whitening operation 603 requires the DoF's of x.sub.0 (k.sub.1,n.sub.sym), M.sub.DoF=K.sub.0, be substantively larger than the number of inner-code stacked signal vectors, N.sub.data=K.sub.1N.sub.sym. If the tapering windows 602 are rectangular, the SINR-revealing metric {circumflex over (γ)}.sub.FLS() can be optionally converted to unbiased SINR estimate

[00030]$\begin{array}{cc}{\tilde{\gamma }}_{\mathrm{FLS}}\left(\ell \right)=\left(1-\frac{K_0-1}{K_1{N}_{\mathrm{sym}}}\right){\hat{\gamma }}_{\mathrm{FLS}}\left(\ell \right)-\frac{K_0}{K_1{N}_{\mathrm{sym}}},& \left(\mathrm{Eq36}\right)\end{array}$

which holds closely if the interference is i.i.d. complex-Gaussian over the outer-code dimension.

[0130] Assuming uniform FOA spacing α(k.sub.FOA)=(k.sub.FOA/K.sub.FOA)f.sub.sym, (Eq28) can be computed using K.sub.0K.sub.1=K.sub.sub N.sub.sym:K.sub.FOA efficient DFT operations. Similarly, assuming uniform TOA spacing τ(n.sub.TOA)=(n.sub.FOA/N.sub.TOA)T.sub.sym/K.sub.0, (Eq29) can be computed using K.sub.0K.sub.FOAL K.sub.1:M.sub.TOA efficient inverse-DFT (IDFT) operations. These operations are both highly regular and parallelizable, allowing their implementation using efficient FPGA or general-purpose GPU (GPGPU) computation modules.

[0131] The FLS TOA-FOA estimates given in (Eq33), and the FLS SINR, are further refined using local search methods in the vicinity of the maximizing FLS TOA-FOA surface value 508. Simple methods for accomplishing include polynomial fit to the surface peak, e.g., using two-dimensional quadratic fit over the nearest neighbors to the maximizing surface grid location. Optionally, parametric search operations that exploit the fully-dispersive form of d.sub.1(k.sub.1,n.sub.sym; τ, α; ) given in (Eq25) multiplied by δ.sub.1(k.sub.1,n.sub.sym;α) in (Eq27), can be used. In one aspect, Newton and Gauss-Newton recursions are defined over the symbol-normalized TOA-FOA vector

[00031]$\begin{array}{cc}{\hat{\upsilon }}_1\left(\ell \right)=\left(\begin{array}{c}{\hat{\tau }}_{T\left(\ell \right)R}{K}_0{f}_{\mathrm{sym}}-\mathrm{floor}\left({\hat{\tau }}_{T\left(\ell \right)R}{K}_0{f}_{\mathrm{sym}}\right)\\ {\hat{\alpha }}_{T\left(\ell \right)R}{T}_{\mathrm{sym}}-\mathrm{round}\left({\hat{\alpha }}_{T\left(\ell \right)R}{T}_{\mathrm{sym}}\right)\end{array}\right).& \left(\mathrm{Eq37}\right)\end{array}$

Defining 2×1 symbol-normalized frequency-time vector

[00032]$\begin{array}{cc}{g}_1\left(k_1,n_{\mathrm{sym}}\right)\stackrel{△}{=}\left(\begin{array}{c}-\left(k_1-\frac{K_1-1}{2}\right)\\ \left(1+\frac{f_1\left(k_1\right)}{f_T}\right)\left(n_{\mathrm{sym}}-\frac{N_{\mathrm{sym}}-1}{2}\right)\end{array}\right),\{\begin{array}{c}0\le k_1<K_1\\ 0\le n_{\mathrm{sym}}<N_{\mathrm{sym}}\end{array},& \left(\mathrm{Eq38}\right)\end{array}$

and K.sub.1N.sub.sym×2 symbol-normalized matrix G.sub.1=[g.sub.1.sup.T(k.sub.1,n.sub.sym)], the recursion is given by

[00033]$\begin{array}{cc}{d}_1\left(\ell \right)=\left[\sqrt{w_{\mathrm{TOA}}\left(k_1\right)w_{\mathrm{FOA}}\left(n_{\mathrm{sym}}\right)}\exp \left\{j({\theta }_{T_1\left(\ell \right)}\left(k_1\right)+2\pi g_1^T\left(k_1,n_{\mathrm{sym}}\right){\hat{\upsilon }}_1\left(\ell \right)\right\}\right],& \left(\mathrm{Eq39}\right)\end{array}$$\begin{array}{cc}{u}_{\mathrm{FLS}}\left(\ell \right)=Q_0^{H}d\left(\ell \right),& \left(\mathrm{Eq40}\right)\end{array}$$\begin{array}{cc}{y}_1\left(\ell \right)=Q_0{u}_{\mathrm{FLS}}\left(\ell \right),& \left(\mathrm{Eq41}\right)\end{array}$$\begin{array}{cc}{g}_{\mathrm{FLS}}\left(\ell \right)=G_1^T\mathrm{Im}\left\{y_1\left(\ell \right)\circ d_1^{*}\left(\ell \right)\right\},& \left(\mathrm{Eq42}\right)\end{array}$$\begin{array}{cc}{H_{\mathrm{FLS}}\left(\ell \right)}_{+}=G_1^T\mathrm{Re}\left\{\left(Q_0{Q}_0^H\right)\circ \left(d_1^{*}\left(\ell \right)d_1^T\left(\ell \right)\right)\right\}{G}_1,& \left(\mathrm{Eq43}\right)\end{array}$$\begin{array}{cc}{H}_{\mathrm{FLS}}\left(\ell \right)={H_{\mathrm{FLS}}\left(\ell \right)}_{+}+G_1^T\mathrm{diag}\left\{\mathrm{Re}\left(y_1\left(\ell \right)\circ d_1^{*}\left(\ell \right)\right)\right\}{G}_1,& \left(\mathrm{Eq44}\right)\end{array}$$\begin{array}{cc}{\hat{\upsilon }}_1\left(\ell \right)\leftarrow \{\begin{array}{cc}{\hat{\upsilon }}_1\left(\ell \right)+\frac{1}{2\pi }{H_{\mathrm{FLS}}^{-1}\left(\ell \right)}_{+}{g}_{\mathrm{FLS}}\left(\ell \right),& \mathrm{Gauss}-\mathrm{Newton}\mathrm{recursion}\\ {\hat{\upsilon }}_1\left(\ell \right)+\frac{1}{2\pi }{H}_{\mathrm{FLS}}^{-1}\left(\ell \right)g_{\mathrm{FLS}}\left(\ell \right),& \mathrm{Newton}\mathrm{recursion}\end{array}& \left(\mathrm{Eq45}\right)\end{array}$$\begin{array}{cc}{\hat{\upsilon }}_1\left(\ell \right)\leftarrow \left(\begin{array}{c}{\left({\hat{\upsilon }}_1\left(\ell \right)\right)}_1-\mathrm{floor}\left\{{\left({\hat{\upsilon }}_1\left(\ell \right)\right)}_1\right\}\\ {\left({\hat{\upsilon }}_1\left(\ell \right)\right)}_2-\mathrm{round}\left\{{\left({\hat{\upsilon }}_1\left(\ell \right)\right)}_2\right\}\end{array}\right).& \left(\mathrm{Eq46}\right)\end{array}$

The final TOA and FOA are then given =({circumflex over (υ)}.sub.1()).sub.1T.sub.sym/K.sub.0 and =({circumflex over (υ)}.sub.1()).sub.2f.sub.sym, and the FLS SINR estimate is given by {circumflex over (γ)}.sub.FLS()=∥u.sub.FLS()∥.sub.2.sup.2/1−∥u.sub.FLS()∥.sub.2.sup.2).

[0132] Equation (Eq33) shows that the FLS TOA-FOA spectrum 613 possesses TOA and FOA ambiguity T.sub.sym/K.sub.0 and f.sub.sym, respectively. In some aspects, this ambiguity is resolved using copy-aided parameter estimation methods 510 that exploit the model of a.sub.0 (τ;) given in (Eq24). In one aspect, the copy-aided ambiguity resolution algorithm is given by

[00034]$\begin{array}{cc}{u}_0\left(n_{\mathrm{zone}},k_{\mathrm{tile}};\ell \right)={\left[\exp \left\{-j{\theta }_{T_0\left(\ell \right)}\left(\left(k_0-k_{\mathrm{tile}}\right)\mathrm{mod}{K}_0\right)\right\}\right]}_{k_0=0}^{K_0-1}\circ {\left[\exp \left\{j2\pi f_0\left(k_0\right)\left({\hat{\tau }}_{T\left(\ell \right)R}+n_{\mathrm{zone}}\frac{T_{\mathrm{sym}}}{K_0}\right)\right\}\right]}_{k_0=0}^{K_0-1},& \left(\mathrm{Eq47}\right)\end{array}$$\begin{array}{cc}{v}_0\left(n_{\mathrm{zone}},k_{\mathrm{tile}};\ell \right)=C_0^H{u}_0\left(n_{\mathrm{zone}},k_{\mathrm{tile}};\ell \right),& \left(\mathrm{Eq48}\right)\end{array}$$\begin{array}{cc}{g}_0\left(n_{\mathrm{zone}},k_{\mathrm{tile}};\ell \right)=u_{\mathrm{FLS}}^H\left(\ell \right)v_0\left(n_{\mathrm{zone}},k_{\mathrm{tile}};\ell \right),& \left(\mathrm{Eq49}\right)\end{array}$$\begin{array}{cc}{S}_{\mathrm{FLS}}\left(n_{\mathrm{zone}},k_{\mathrm{tile}};\ell \right)={.Math.g_0\left(n_{\mathrm{zone}},k_{\mathrm{tile}};\ell \right).Math.}^2/{.Math.v_0\left(n_{\mathrm{zone}},k_{\mathrm{tile}};\ell \right).Math.}_2^2,& \left(\mathrm{Eq50}\right)\end{array}$$\begin{array}{cc}\left\{{\hat{n}}_{\mathrm{zone}}\left(\ell \right),{\hat{k}}_{\mathrm{tile}}\left(\ell \right)\right\}=\arg \underset{\begin{array}{c}{n}_{\mathrm{zone}}=0,\rightleftharpoons ,K_0-1,\\ k_{\mathrm{tile}}=-K_{\mathrm{tile}},\rightleftharpoons ,+K_{\mathrm{tile}}\end{array}}{\max }{S}_{\mathrm{FLS}}\left(n_{\mathrm{zone}},k_{\mathrm{tile}};\ell \right),& \left(\mathrm{Eq51}\right)\end{array}$

where u.sub.FLS() is given by (Eq40), and where {circumflex over (n)}.sub.zone() and {circumflex over (k)}.sub.tile() are the TOA zone and FOA tile containing the beacon detection, respectively. The full TOA and FOA geo-observables are then given by

[00035]${\hat{\tau }}_{T\left(\ell \right)R}\leftarrow {\hat{\tau }}_{T\left(\ell \right)R}+{\hat{n}}_{\mathrm{zone}}\left(\ell \right)\frac{T_{\mathrm{sym}}}{K_0}\mathrm{and}{\hat{\alpha }}_{T\left(\ell \right)R}\leftarrow {\hat{\alpha }}_R+{\hat{\alpha }}_{T\left(\ell \right)R}+{\hat{k}}_{\mathrm{tile}}\left(\ell \right)f_{\mathrm{sym}},$

respectively. If {circumflex over (k)}.sub.tile()≠0 for any detected beacon, then the FOA centroid 612 {circumflex over (α)}.sub.R is optionally recomputed, e.g., using weighted estimate

[00036]${\hat{\alpha }}_R\leftarrow \underset{\ell \in {\ell }_{\mathrm{detect}}}{.Math.}{\hat{\gamma }}_{\mathrm{FLS}}\left(\ell \right){\hat{\alpha }}_{T\left(\ell \right)R}/\underset{\ell \in {\ell }_{\mathrm{detect}}}{.Math.}{\hat{\gamma }}_{\mathrm{FLS}}\left(\ell \right),$

and the channelized snapshot 553 is regenerated 505 and the subsequent FLS geo-observable estimation operations shown in FIG. 12 are repeated with the new FOA centroid 612 estimate.

[0133] The copy-aided ambiguity estimator also provides complex gain estimate

[00037]$\begin{array}{cc}{\hat{a}}_{\mathrm{sub}}\left(\ell \right)=\frac{g_0\left({\hat{n}}_{\mathrm{zone}}\left(\ell \right),{\hat{k}}_{\mathrm{tile}}\left(\ell \right);\ell \right)}{\begin{array}{c}{.Math.g_0\left({\hat{n}}_{\mathrm{zone}}\left(\ell \right),{\hat{k}}_{\mathrm{tile}}\left(\ell \right);\ell \right).Math.}^2+\\ \left(1-{.Math.u_{\mathrm{FLS}}\left(\ell \right).Math.}_2^2\right){.Math.v_0\left({\hat{n}}_{\mathrm{zone}}\left(\ell \right),{\hat{k}}_{\mathrm{tile}};\ell \right).Math.}_2^2\end{array}},& \left(\mathrm{Eq52}\right)\end{array}$

which can be used to compute the beacon channelizer output power and phase offset. These parameters provide key inputs for channel calibration operations, e.g., to determine transmit and receive carrier phase, and true channel pathloss, for subsequent network calibration operations.

[0134] Table 3 lists FLS surface generation parameters usable in the exemplary UTM and IIoT scenarios described here. The degrees of freedom are large enough to separate all of the beacons in the users' FoV's for each scenario, which a roughly factor-of-two excess to account for multipath reflections. In each case, the number of data entries is a large enough multiple of the despreader DoF's to yield a stable QRD and FLS estimate.

TABLE-US-00003 TABLE 3 Exemplary FLS Surface Generation Parameters Parameter LMS Scenario 2.4 GHz Ch. 13 Scenario K.sub.0 (M.sub.DoF) 32 80 K.sub.1N.sub.sym (N.sub.data) 555 2,200 K.sub.FOA 128 bins 128 bins M.sub.TOA 32 lags 128 lags

[0135] Table 4 summarizes complexity of the channelization and FLS surface generation operations for the two scenarios, assuming one real multiply-and-add per operation, and assuming the whitening operation 603 is performed using a QRD instantiated using Modified Gram-Schmidt Orthogonalization (MGSO). Of these operations, the QRD uses than 23% of the total operations in each scenario. The complexity is well within the capabilities of modern DSP gear, even for the 2.4 GHz Ch. 13 Scenario. Moreover, the channelization and DFT/IDFT operations are easily implemented in FPGA (e.g., Xilinx 7K325T or higher devices) or using general-purpose GPU's (GPGPU's).

TABLE-US-00004 TABLE 4 Channelization, FLS TOA-FOA Surface Generation Complexity Operation LMS Scenario 2.4 GHz Ch. 13 Scenario Channelization 7.06 Mops 61.9 Mops FLS whitening 2.27 Mops 56.8 Mops FLS surface 6.57 Mops 130.2 Mops Total FLS 15.9 Mops 248.9 Mops

[0136] Table 5 summarizes memory requirements of the FLS surface generation procedure for the two scenarios. The memory requirements assuming 64-bit data precision at each stage of processing, and in-place QRD operations. The channelization operation imposes the bulk of memory requirements for each scenario, and is not particularly onerous in any case. The channelization memory is also well within capability of modern DSP, FPGA (e.g., Xilinx 7K325T or higher devices), or GPGPU's.

TABLE-US-00005 TABLE 5 Channelization, FLS TOA-FOA Surface Generation Memory Operation LMS Scenario 2.4 GHz Ch. 13 Scenario Channelization 2 MB 15.54 MB FLS whitening 21 KB 119 KB FLS surface 11 KB 64 KB Total FLS 2.03 MB 15.72 MB

[0137] FIG. 16A FIG. 16B, FIG. 17A and FIG. 17B show the FLS TOA-FOA surface 613 obtaining for the LMS Scenario and 2.4 GHz Ch. 13 Scenario, respectively, and for the surface instantiation parameters given in Table 3. Beacon 10 (out of 32 beacons searched, and 18 beacons detected, i.e., every network node 401 in the FoV of User 1 400) is chosen for the LMS Scenario, and Beacon 37 (out of 50 beacons searched, and 29 beacons detected, i.e., every network node 401 in the FoV of User 1 400) is chosen for the 2.4 GHz Ch. 13 Scenario, because they each provide a median peak value over all the beacons detected during their respective TOA-FOA searches.

[0138] FIG. 16A shows the exemplary LMS Scenario FLS TOA-FOA surface 613 over its full range, generated using the simplified outer-code signal model given in (Eq25). The FLS TOA-FOA surface 613 also shows detection of TOA zone 2, using the copy-aided TOA-FOA ambiguity resolution method 510 described in (Eq47)-(Eq51). The FOA clutter spectrum 610 is clearly visible in the full TOA-FOA surface 613, as well as a single TOA-FOA peak obtained at a 14.8 dB estimated SINR, or 17.1 dB above the background clutter. No other peaks associated with any other network node 401 are visible, despite the fact that 18 network nodes 401 are in the user's 400 FoV.

[0139] FIG. 16B shows results of a 500×500 parametric local search 508 centered on the peak of the full FLS TOA-FOA surface 613 shown in FIG. 16A, and generated using the extended outer-code signal model given by (Eq25), multiplied by dispersive term (Eq27). The benefit of the extended outer-code signal model is evident in the maximal peak value of 50.4 dB—a 35.6 dB increase over the full TOA-FOA surface 613, and 52.8 dB above the FOA clutter spectrum 610. The TOA-FOA error for this case is (7.8 μs, −16.7 mHz), comparable to the (51.6 μs, 10.6 mHz) 95.sup.th percentile CRB for this case, showing that the method is adhering well to expected performance.

[0140] FIG. 17A shows the exemplary 2.4 GHz Ch. 13 FLS TOA-FOA surface 613 over its full range, generated using the simplified outer-code signal model given in (Eq25). The FLS TOA-FOA surface 613 also shows detection of TOA zone 59, using the copy-aided TOA-FOA ambiguity resolution method 510 described in (Eq47)-(Eq51). The FOA clutter spectrum 610 is clearly visible in this surface, as well as a single TOA-FOA peak obtained at an 11.4 dB estimated SINR, 19.4 dB above the background clutter. No other peaks are visible, despite the fact that 29 network nodes 401 are in the user's 400 FoV.

[0141] FIG. 17B shows a 50×50 parametric local search 508, centered on the peak of the full FLS TOA-FOA surface 613 shown in FIG. 17A, and generated using the extended outer-code signal model given by (Eq25), multiplied by dispersive term (Eq27). The benefit of the extended outer-code signal model is also evident in the maximal peak value of 25 dB—a 13.6 dB increase over the full FLS TOA-FOA surface 613, and 33 dB above the FOA clutter spectrum 610. The TOA-FOA error for this case is (−7.9 μs, 16.1 mHz), comparable to the (39.6 μs, 82.2 mHz) 95.sup.th percentile CRB for this case, showing that the method is also adhering well to expected performance in this Scenario.

[0142] Other aspects of the disclosure employ similar operations using alternate stacking methods. These include outer-code stacking, which transforms x.sub.sub(k.sub.sub,n.sub.sym) into K.sub.1×1 outer-code stacked vector x.sub.1(k.sub.0,n.sub.sym)=[x.sub.sub(K.sub.0k.sub.1+k.sub.0,n.sub.sym)].sub.k.sub.1.sub.=0.sup.K.sup.1.sup.−1, and symbol stacking, which transforms x.sub.sub(k.sub.sub,n.sub.sym) into N.sub.sym×1 symbol-stacked vector x.sub.sym(k.sub.sub,n.sub.sym)=[x.sub.sub(k.sub.sub,n.sub.sym)].sub.n.sub.sym.sub.=0.sup.N.sup.sym.sup.−1. Each stacking operation has advantages and disadvantages. Inner-code stacking results in a TOA estimate ≈mod(T.sub.sym/K.sub.0), i.e., with range aliased between 0 and T.sub.sym/K.sub.0, but with high precision within that range. For this reason it is referred to here as the fine least-squares (FLS) estimator. Conversely, outer-code stacking yields TOA estimate ≈mod T.sub.sym, i.e., with range aliased to full range between 0 and T.sub.sym, but with precision that is a factor of K.sub.0 coarser. It is referred to here as the coarse least-squares (CLS) estimator. Both estimators provide full range and precision in FOA. The symbol-stacked estimator provides full range and precision in TOA, but no estimate of FOA estimate. In all cases, these issues are resolvable using copy-aided post-processing methods 508.

Geo-Observable Based Positioning/Timing Procedure

[0143] FIG. 18 shows the geo-observable based positioning and timing procedure disclosed here, applicable to at least the network shown in FIG. 1, and describes the common operations used in each stage. The procedure optimizes maximum-likelihood (ML) objective function

[00038]$\begin{array}{cc}{F}_{\mathrm{MLE}}\left(p_R,v_R,{\tau }_R,{\alpha }_R\right)=\frac{1}{2}\left(\frac{1}{{\sigma }_{\mathrm{TOA}}^2}{F}_{\mathrm{TOA}}\left(p_R,{\tau }_R\right)+\frac{1}{{\sigma }_{\mathrm{FOA}}^2}{F}_{\mathrm{FOA}}\left(p_R,v_R,{\alpha }_R\right)\right),& \left(\mathrm{Eq53}\right)\end{array}$

where F.sub.TOA (p.sub.R, τ.sub.R) and F.sub.FOA (p.sub.R,v.sub.R,α.sub.R) are the TOA-only and FOA-only ML estimators, respectively,

[00039]$\begin{array}{cc}{F}_{\mathrm{TOA}}\left(p_R,{\tau }_R\right)=\underset{\ell \in {ℒ}_{\mathrm{detect}}}{.Math.}{\hat{\gamma }}_{\left(\right)}\left(\ell \right){\left({\hat{\tau }}_{T\left(\ell \right)R}-\left({\tau }_R+\frac{1}{c}{.Math.p_{T\left(\ell \right)R}.Math.}_2\right)\right)}^2,& \left(\mathrm{Eq54}\right)\end{array}$$\begin{array}{cc}{F}_{\mathrm{FOA}}\left(p_R,v_R,{\alpha }_R\right)\approx \underset{\ell \in {ℒ}_{\mathrm{detect}}}{.Math.}{\hat{\gamma }}_{\left(\right)}\left(\ell \right){\left({\hat{\alpha }}_{T\left(\ell \right)R}-\left({\alpha }_R+\frac{1}{{\lambda }_T}{u}_{T\left(\ell \right)R}^T{v}_R\right)\right)}^2,& \left(\mathrm{Eq55}\right)\end{array}$

and where λ.sub.T=f.sub.T/c is the nominal signal-in-space wavelength, σ.sub.TOA.sup.2 and σ.sub.FOA.sup.2 are given in (Eq34) and (Eq35), respectively. In the aspect described here, the networks nodes 401 are assumed to be fixed and have known positions, such that network node 401 position (t.sub.UTC) ≡ and velocity (t.sub.UTC)≡0. If the network nodes' 401 system clocks 214 are synchronized to UTC, then the user 400 system clock 214 and LO offsets in the reception operations 300 are given by τ.sub.R=t.sub.R.sup.ref−t.sub.UTC(t.sub.R.sup.ref)=t.sub.R.sup.ref−t.sub.UTC.sup.ref and α.sub.R=−f.sub.Tε.sub.R/(1+ε.sub.R), respectively, allowing the UTC time at known receive reference time t.sub.R.sup.ref and the clock rate offset ε.sub.R to be derived from the observed estimates, and the user 400 position and velocity being estimated by the method are given by p.sub.R=p.sub.R(t.sub.UTC.sup.ref) and v.sub.R=v.sub.R(t.sub.UTC.sup.ref), respectively.

[0144] In some aspects, the network nodes' 401 system clocks 214 are not fully synchronized to UTC, but network node 401 timing offsets =−t.sub.UTC and rate offsets from UTC have been estimated, e.g., using network calibration procedures. In this case, the detected TOA's and FOA's are adjusted to compensate for these offsets. In one aspect, this is performed by setting

←−(−(1+){circumflex over (t)}.sub.UTC()),  (Eq56)

←+f.sub.T,  (Eq57)

prior to computation of the ML estimator. It should be noted that (Eq57) is not exact, as it fails to include division of the second term by 1+ε.sub.R as shown in (Eq14). However, this effect is minor for user 400 system clocks 214 with <20 ppm rate offset, and can be removed in subsequent refinements.

[0145] Estimating and concentrating τ.sub.R and α.sub.R estimates out of (Eq54)-(Eq55) yields

[00040]$\begin{array}{cc}{\hat{\tau }}_R{.Math.}_{\mathrm{ML}}=\underset{\ell \in {ℒ}_{\mathrm{detect}}}{.Math.}\left(\frac{{\hat{\gamma }}_{\left(\right)}\left(\ell \right)}{\underset{{\ell }^{\prime }\in {ℒ}_{\mathrm{detect}}}{.Math.}{\hat{\gamma }}_{\left(\right)}\left({\ell }^{\prime }\right)}\right)\left({\hat{\tau }}_{T\left(\ell \right)R}-\frac{1}{c}{.Math.p_{T\left(\ell \right)R}.Math.}_2\right)={.Math.{\hat{\tau }}_{T\left(\ell \right)R}-\frac{1}{c}{.Math.p_{T\left(\ell \right)R}.Math.}_2.Math.}_{{ℒ}_{\mathrm{detect}}},& \left(\mathrm{Eq58}\right)\end{array}$$\begin{array}{cc}.Math.F_{\mathrm{TOA}}\left(p_R\right)\leftarrow {.Math.{\left({\hat{\tau }}_{T\left(\ell \right)R}-\frac{1}{c}{.Math.p_{T\left(\ell \right)R}.Math.}_2\right)}^2.Math.}_{{ℒ}_{\mathrm{detect}}}-{.Math.{\hat{\tau }}_{T\left(\ell \right)R}-\frac{1}{c}{.Math.p_{T\left(\ell \right)R}.Math.}_2.Math.}_{{ℒ}_{\mathrm{detect}}}^2={.Math.{\left(\left({\hat{\tau }}_{T\left(\ell \right)R}-{.Math.{\hat{\tau }}_{T\left(\ell \right)R}.Math.}_{{ℒ}_{\mathrm{detect}}}\right)-\frac{1}{c}\left({.Math.p_{T\left(\ell \right)R}.Math.}_2-{.Math.{.Math.p_{T\left(\ell \right)R}.Math.}_2.Math.}_{{ℒ}_{\mathrm{detect}}}\right)\right)}^2.Math.}_{{ℒ}_{\mathrm{detect}}}.& \left(\mathrm{Eq59}\right)\end{array}$$\begin{array}{cc}{\hat{\alpha }}_R{.Math.}_{\mathrm{ML}}={.Math.{\hat{\alpha }}_{T\left(\ell \right)R}-\frac{1}{{\lambda }_T}{u}_{T\left(\ell \right)R}^T{v}_{T\left(\ell \right)R}.Math.}_{{ℒ}_{\mathrm{detect}}},{\lambda }_T=\frac{c}{f_T},& \left(\mathrm{Eq60}\right)\end{array}$$\begin{array}{cc}.Math.F_{\mathrm{FOA}}\left(p_R,v_R\right)\leftarrow .Math.{\left({\hat{\alpha }}_{T\left(\ell \right)R}-\frac{1}{{\lambda }_T}{u}_{T\left(\ell \right)R}^T{v}_R\right)}^2.Math.-{.Math.{\hat{\alpha }}_{T\left(\ell \right)R}-\frac{1}{{\lambda }_T}{u}_{T\left(\ell \right)R}^T{v}_R.Math.}_{{ℒ}_{\mathrm{detect}}}^2=.Math.{\left(\left({\hat{\alpha }}_{T\left(\ell \right)R}-{.Math.{\hat{\alpha }}_{T\left(\ell \right)R}.Math.}_{{ℒ}_{\mathrm{detect}}}\right)-\frac{1}{{\lambda }_T}{\left(u_{T\left(\ell \right)R}-.Math.u_{T\left(\ell \right)R}-{.Math.u_{T\left(\ell \right)R}.Math.}_{{ℒ}_{\mathrm{detect}}}\right)}^T{v}_R\right)}^2.Math..& \left(\mathrm{Eq61}\right)\end{array}$