Julia and the GPS payload onboard Waratah Seed-1 satellite

Community
1

Hi Julia community,

I’m Ignatius Rivaldi, research engineer working at UNSW ACSER: Australian Centre for Space Engineering Research (ACSER)

We flew our GPS receiver payload, Harry v3 on Waratah Seed-1 6U cubesat developed by our collaborator, ARC Training Centre for CubeSats, UAVs, and their Applications (CUAVA)

After several public announcements and a presentation at ASRC 2024 , this is our first data from our Kea GPS receiver navigating in space:

|794px;x435px;
|794px;x435px;1920×1052 926 KB

(The TLE was propagated by SatelliteToolbox.jl and plotted by Makie.jl)

LinkedIn post by our PI, Andrew Dempster on LinkedIn: To you, this might look like an ordinary image of a GNSS receiver… | 23 comments

Why Julia

Part of the process of hardware integration between our GPS payload and WS-1 bus is verifying whether the RF chain works before we launch it into space. In orbit, our FPGA-based Kea GPS receiver will process the GPS signals into navigation solutions. But, if something breaks, it will be hard to debug.

Our payload can be configured to dump raw intermediate frequency (IF) samples from its RF frontends at data rate of about 2 MB/s . We can then run a software receiver on the IF samples to verify that our RF frontends can successfully receive GPS signal. A good guide of this can be found in Daniel Estevez’s post: Timing SDR recordings with GPS – Daniel Estévez

There are a lot of open source GPS software receivers out there:

This is based on GNU Radio framework , and this is the most featureful open source GNSS software receiver out there. It can process GPS L1,L2,L5, Galileo, Beidou, GLONASS including Galileo navigation message authentication,
at real time on embedded hardware such as Raspberry Pi. The issue is that to support those functions and as it is using GNU Radio, it is written in some extremely involved C++: Fundamentals - GNSS-SDR

This is written in Matlab, have similar functions as GNSS-SDR. Because it is written in Matlab, it is much easier to modify and extend compared to GNSS-SDR.

Unfortunately, Matlab doesn’t encourage the user to write performant code, from hiding how their JIT compiler worked, to having multithreading support as a paid $1800 DLC. The worst sin is that in 2024, Matlab doesn’t have an official function of reading interleaved complex binary file, leading to Matlab users using some very inefficient and error-prone method of reading complex data or third-party solution that implements complex file read in C.

All of these papercuts leads to this software receiver taking 30 minutes! to process 40 second of IF data, while GNSS-SDR runs in only 10 seconds

We have a very fast software receiver that’s hard to extend, and a very slow software receiver that’s easy to extend. Sounds like a two-language problem to me

Then I found this from @zsoerenm :

This is perfect for my application,as I can:

  • Read any complex file is just a mmap("complex.iq", Complex{Float32}) away, no matter what the real datatype is, from UInt8 to Float64, with it just working by the power of parametric polymorphism and multiple dispatch
  • see the assembly output of Julia compiler by @code_native,
  • benchmark by BenchmarkTools.jl,
  • profile by VS code @profview ,
  • use Tracy.jl ,the profiler that game devs uses
  • use LIKWID.jl to measure FLOPS
  • use Bumper.jl to reduce memory allocations and GC
  • Threads.@threads for is free!

Now we have a GPS software receiver that’s as easy to modify as the Matlab receiver, but nearly as fast as GNSS-SDR

Julia running in space!

Julia can also run on Raspberry Pi CM4, the processor I used on our GPS payload computer. Unfortunately, precompilation delay means that its not practical to run the GPS signal processing script in space - every watt matters in space.

In space, we also capture IF data for us to do postprocessing in ground. But the IF data doesn’t fit into a single 1-2 minute satellite downlink pass, so I wrote a small Julia script to chunk the IF file on the payload computer into multiple smaller chunks, which can be downloaded over multiple passes. Which is then recombined at the ground, then postprocessed using JuliaGNSS

|  WS-1 GPS payload |   WS-1 cubesat bus       |   Groundstation
Kea -> Julia script ->  S-band radio to ground   ->  JuliaGNSS

This is the output of one of our IF captures that was chunked by Julia script in space, and postprocessed by JuliaGNSS at ground:

|292px;x194px;
|292px;x194px;593×394 38.7 KB

JuliaGNSS can receive and decode 3 GPS satellites:

|415px;x97px;
|415px;x97px;516×199 7.06 KB

This is the first time I know of that Julia is running in space

89 Likes
2

It’s great to see that JuliaGNSS is actually being used :wink:
It motivates me to continue my work.

Now we have a GPS software receiver that’s as easy to modify as the Matlab receiver, but nearly as fast as GNSS-SDR

I have more performance improvements in the pipe :slight_smile:

19 Likes
3

It’s really excellent!! Note that if you want to connect your processing more directly with software-defined radios, you can take a look at the SDR driver we developed on JuliaTelecom to directly output the IQ data to your processing unit. But it’s truly exciting to see what can be achieved with Julia-based signal processing!

7 Likes
4

There are more reasons that makes Julia one of the very promising language and ecosystem for SDR signal processing:

There are PLL algorithms that promises to be ‘better’ than conventional second order PLL loop:

But GNU Radio is still stuck in second order PLL and symbol synchronization loop. What gives?

One of the possible reason is that GNU Radio 3.x flowgraph doesn’t support feedback: Why can’t we do loops?

A lot of users come into the project looking to experiment, including building known algorithms like phase-locked loops (PLLs). GNU Radio comes with each of the blocks to build your own PLL in a flowgraph, like a multiplier, lowpass filter, and VCO. But you try to put it together and GNU Radio tells you it won’t work. That the flow graphs don’t allow you to do loops. Why? PLLs are fundamental to radio and signal processing! How can it be that I can’t build loops into a flow graph?

So a user that wants to implement their exotic PLLs to evaluate whether it meets their needs can’t use the user friendly GNU Radio Companion flowgraph interface:

image
image1668×1061 214 KB

and they need to drop to some involved C++:

It turns out that GNU Radio consists of 4 languages:

Its hard to move between dragging virtual cables between signal processing blocks and writing C++ .

While Julia might be able to replace the 4 different languages with one using metaprogramming and code generation.

It seems that they simplify this situation in GNU Radio 4.0 using modern C++. I haven’t played around with it though.

Another trend is because of AI hype, everyone is making their own AI accelerator chip. Those are screaming for us to run SDR workloads on that, all AI accelerator chip do is accelerate the matmuls. They don’t care whether its part of deep neural net or signal processing flowgraph.
Tenstorrent have very interesting architecture:

image
image2334×1314 510 KB

What we need is a good framework to convert our SDR flowgraphs into something that can run on Tensix core array. GNU Radio can’t do that because they connect precompiled C++ blocks, while Julia GPUCompiler.jl might be able to compile our flowgraph description into Tensix RISC-V assembly.

6 Likes
5

Could the first OP post be modified to become a Julia-lang blog post?

I wish Julia-lang blog could also post successful use cases of Julia like this one.

3 Likes
6

Oh that is a great idea. I am pinging @rayegun !

1 Like
7

Hi @minetest2048 !

This is awesome! I am so glad to see more people using Julia in the space community :slight_smile:

One of my projects now is to prove that we can use Julia for RT applications to write an entire attitude and orbit control subsystem software of a satellite in Julia. I am still having problems for RT deadlines 5% of the time, but when we fix that, I will apply for funding a CubeSat mission to test the approach. In this case, I will mention your mission that now has Julia running in space. It will be some much better for the proposal!

Thank you very much for sharing!

7 Likes
8

I should get access to a tenstorrent board in the next few months and I want to experiment with julia on it. But I was planning to go by the MLIR route: GitHub - tenstorrent/tt-mlir: Tenstorrent MLIR compiler

5 Likes
9

Thanks @Ronis_BR for your SatelliteToolbox.jl, I used it a lot for mission planning

Some wishlist for SatelliteToolbox:

  • I wish we have support for converting ECI to RIC frame, all the existing literature uses RIC frame to measure error between GPS position fixes and propagated TLE, as it is easier to interpret than ECEF error
  • Typed data input and output so that I can just compare the output of TLE propagation with GPS fixes from our receiver without messing with Julian days so I can convert ECI to ECEF or trying to understand which 100 reference frames I should use so my TLE works correctly (turns out its PEF)

Ideally it works like this:

orbp = Propagators.init(Val(:SGP4), iss_tle)
ecis = Propagators.propagate!.(orbp, 0:1:86400)
ecefs = to_ecef(ecis) # It handles all the conversions without me needing to construct the transform matrix

With humans landing on the moon in hopefully 5 years, things will become more interesting. We will have lunar coordinate system, lunar UTC, lunar GPS. We also going to try to receive Earth GPS on the Moon and try to navigate with that using NaviMoon receiver

For your realtime Julia, this CppCon talk just dropped yesterday about Clang RealtimeSanitizer and performance constraints for C++:

Previous thread on this: Can we get RealtimeSanitizer in Julia?

If we take inspiration from this, and extend AllocCheck.jl to not just detect memory allocation, but also other non-RT functions such as gc_collect and taking locks/mutexes, we can go closer to realtime Julia.

The other option that works today, but its very brutal is compiling your code to a shared library using StaticCompiler.jl. You won’t trigger GC if there’s no GC in the first place, same with memory allocation.

I tried to write a GNU Radio block with this and kind of gave up. Its very easy to accidentally depend on libjulia and have missing symbols errors during linking. But the worst thing is that StaticCompiler doesn’t currently emit debug symbols. So when you segfault, which is very easy to do because freestanding Julia doesn’t currently have memory safety features unlike Rust, both gdb and Valgrind won’t help you find where it is and debug it

2 Likes
10

CoordRefSystems.jl might help with this wishlist item if I understood the challenge correctly.

1 Like
11

Hi @minetest2048 !

Awesome!

I can easiy implement this. Do you mind open an issue so that I do not forget?

We can do this using a layer on top of the functions we have now. However, we must think very careful because, depending on the use case, you might have some problems. There are a lot of ECEF reference frames. You selected PEF from IAU-76 convention, but another application might require the ITRF or another one from the IAU-2010 convertion.

That’s why I decided to let the user select the specific frame.

I fully agree!

12

Ah! By the way, one of my projects if to extend static arrays so that we can annotate the reference system in which the vector is represented. It will help a lot, since, in this case, we can convert by only specifying the target system. However, it will require a huge rewrite of a lot of packages in the SatelliteToolbox.jl ecosystem.

13

That is precisely what CoordRefSystems.jl does, except we don’t implement the StaticArrays.jl interface. The AbstractArray interface is mostly useful with Cartesian frames where LinearAlgebra plays a major role.

Also handled in CoordRefSystems.jl? This is the datum type parameter that we store in the type at compile time.

14

Hi @juliohm !

Yes, I know! However, last time I checked, CoordRefSystem.jl could not fit our applications here. We need something more “low-level” when representing reference system transformations.

For example, I think it is not possible to obtain the quaternion, or DCM, that rotates the True-of-Date frame at some epoch to the True-Equator, Mean-Equinox frame in another epoch.

SatelliteToolboxTransformations.jl together with ReferenceFrameRotations.jl were designed since the beginning to take into account the features we need to simulate satellite attitude dynamics and also perform tasks related to mission analysis.

Currently, we obtain the rotation entities (DCM, Quaternion, etc) and the user must track the vector representation. What I want to do is the add a layer so that the reference in which each vector is represented is annotated. This modification will not make everything automatic, but it will decrease the probability of errors.

15

Hi @Ronis_BR !

Appreciate it if you could elaborate on this feature request. We do handle dynamic frames and epochs in CoordRefSystems.jl convert methods.

16

Hi @juliohm!

What we need is basically all the transformations described here:

And we need to obtain the transformed vectors and the entities as well (DCM and Quaternions).

17

What’s currently missing from CoordRefSystems.jl representation is time, and time reference frame

Currently I’m working with several time frames/representations:

  • PositionVelocityTime.jl uses AstroTime.TAI from navigation solution
  • Our Kea GPS receiver outputs GPS time as (week_number, tow) tuple in the NMEA message
  • SatelliteToolbox.jl uses Julian day UTC in unitless Float64
  • We print using DateTime UTC

If we’re off by several milliseconds, at 7 km/s my satellite will be off by several hundred meters as shown here:

@zsoerenm please fix your sign!

GPS broadcasts conversion factors between GPS time and UTC, but for my use case its negligible.

The thing with ECI and ECEF conversion is that the conversion DCM changes every second, unlike geodetic datum conversion where the agency releases new datum every year or so. At least that’s what I understand from SatelliteToolbox ISS tutorial

What’s currently missing is a StateVector struct, parameterized on time type and reference frame type for position and its derivatives:

struct StateVector{T,P,V,reference_frame}
     time::T
     position::P 
     velocity::V 
    #higher order derivatives if necessary 
end

This is more tricky then I first thought…

1 Like
18

To add a blog post, just create a pull request:

4 Likes
19

Yes! You are completely right.

We really need to annotate the reference frame, but we already have a state vector representation that you can use to convert between ECI and ECEF. It also transforms the velocity as it was measured by an observer in the ECEF reference frame:

julia> using SatelliteToolbox

julia> f = create_tle_fetcher(CelestrakTleFetcher)
CelestrakTleFetcher("https://celestrak.org/NORAD/elements/gp.php")

julia> tles = fetch_tles(f; satellite_name = "ISS (ZARYA)")
[ Info: Fetch TLEs from Celestrak using satellite name: "ISS (ZARYA)" ...
1-element Vector{TLE}:
 TLE: ISS (ZARYA) (Epoch = 2024-12-25T11:54:31.088)

julia> iss_tle = tles[1]
TLE:
                      Name : ISS (ZARYA)
          Satellite number : 25544
  International designator : 98067A
        Epoch (Year / Day) : 24 / 360.49619315 (2024-12-25T11:54:31.088)
        Element set number : 999
              Eccentricity :   0.00056820
               Inclination :  51.63780000 deg
                      RAAN :  87.51690000 deg
       Argument of perigee :   4.67930000 deg
              Mean anomaly : 355.42490000 deg
           Mean motion (n) :  15.50207956 revs / day
         Revolution number : 48822
                        B* :   0.00042581 1 / er
                     ṅ / 2 :   0.00024109 rev / day²
                     n̈ / 6 :            0 rev / day³

julia> orbp = Propagators.init(Val(:SGP4), iss_tle)
OrbitPropagatorSgp4{Float64, Float64}:
   Propagator name : SGP4 Orbit Propagator
  Propagator epoch : 2024-12-25T11:54:31.088
  Last propagation : 2024-12-25T11:54:31.088

julia> r_teme, v_teme = Propagators.propagate!(orbp, 0)
([294198.57439655927, 6.7847714223184595e6, 32.68505834064143], [-4751.8787639694765, 199.06599097725888, 6013.0710696574815])

julia> sv_teme = OrbitStateVector(Propagators.epoch(orbp), r_teme, v_teme)
OrbitStateVector{Float64, Float64}:
  epoch : 2.46067e6 (2024-12-25T11:54:31.088)
      r : [294.199, 6784.77, 0.0326851]  km
      v : [-4.75188, 0.199066, 6.01307]  km/s

julia> sv_eci_to_ecef(sv_teme, TEME(), PEF(), Propagators.epoch(orbp))
OrbitStateVector{Float64, Float64}:
  epoch : 2.46067e6 (2024-12-25T11:54:31.088)
      r : [-6758.72, 662.88, 0.0326851]   km
      v : [-0.408953, -4.24116, 6.01307]  km/s

The last output is the position and velocity vector of the ISS as measured by an observer in the PEF reference frame.

You can also add the acceleration to the state vector if you want and the sv_eci_to_ecef will convert it for you.

3 Likes
20

OrbitStateVector is exactly what I wanted! If only we can get Propagators.propagate! to spit out OrbitStateVector automatically instead of (position,velocity) tuple