diff --git a/src/orbit/tle.jl b/src/orbit/tle.jl index ef4a6c18..ff163f4f 100644 --- a/src/orbit/tle.jl +++ b/src/orbit/tle.jl @@ -21,7 +21,7 @@ # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # ==# -export read_tle +export read_tle, read_tle_from_string """ ### macro parse_value(T, str, line_num) @@ -93,7 +93,7 @@ function compute_checksum(str::AbstractString) end """ - function read_tle(tle_filename::String) + @inline function read_tle(tle_filename::String) Read the TLEs in the file `tle_filename`. @@ -108,10 +108,52 @@ Read the TLEs in the file `tle_filename`. An array with all the TLEs that were parsed. """ -function read_tle(tle_filename::String, verify_checksum::Bool = true) +@inline function read_tle(tle_filename::String, verify_checksum::Bool = true) # Open the file in read mode. file = open(tle_filename, "r") + _parse_tle(file, verify_checksum) +end + +""" + @inline function read_tle_from_string(tles::String, verify_checksum::Bool = true) + @inline function read_tle_from_string(tle_l1::String, tle_l2::String, verify_checksum::Bool = false) + +Parse a set of TLEs in the string `tles` or one TLE with first line `tle_l1` and +second line `tle_l2`. + +# Args + +* `tles`: TLEs that will be parsed. +* `verify_checksum`: (OPTIONAL) If false, then the checksum will not be verified + (**Default** = true). + +* `tle_l1`: First line of the TLE. +* `tle_l2`: Second line of the TLE. + +# Returns + +An array with all the TLEs that were parsed. + +""" +@inline function read_tle_from_string(tles::String, + verify_checksum::Bool = true) + # Convert the string to an IOBuffer and call the function to parse it. + _parse_tle(IOBuffer(tles), verify_checksum) +end + +@inline function read_tle_from_string(tle_l1::String, + tle_l2::String, + verify_checksum::Bool = false) + # Assemble the TLE into one string and call the function to parse it. + read_tle_from_string(tle_l1 * '\n' * tle_l2, verify_checksum) +end + +################################################################################ +# Private Functions +################################################################################ + +function _parse_tle(io::IO, verify_checksum::Bool = true) # State machine to read the TLE. It has three possible states: # # :name -> Satellite name. @@ -158,11 +200,11 @@ function read_tle(tle_filename::String, verify_checksum::Bool = true) line = nothing - while !eof(file) + while !eof(io) # Read the current line, strip white spaces, and skip if it is blank or # is a comment. if !skip_line_read - line = strip(readline(file)) + line = strip(readline(io)) line_num += 1 (isempty(line)) && continue (line[1] == '#') && continue diff --git a/test/orbit/propagators/sgp4.jl b/test/orbit/propagators/sgp4.jl index e826b521..aa3a70c2 100644 --- a/test/orbit/propagators/sgp4.jl +++ b/test/orbit/propagators/sgp4.jl @@ -83,5 +83,35 @@ @test st_sgp4_result[7] ≈ SGP4_results[k,7] atol=1e-9 end end + + # Test using the function `read_tle_from_string`. + tle = read_tle_from_string( + "1 23599U 95029B 06171.76535463 .00085586 12891-6 12956-2 0 2905", + "2 23599 6.9327 0.2849 5782022 274.4436 25.2425 4.47796565123555") + + filename = "./sgp4_tests/aiaa-2006-6753/sgp4_tle_23599_result.txt" + SGP4_results = readdlm(filename; comments=true) + + # Initialize the orbit propagator. + orbp = init_orbit_propagator(Val{:sgp4}, tle[1], sgp4_gc_wgs72) + + # Propagate the orbit. + t = SGP4_results[:,1]*60 + (orbm, r_TEME, v_TEME) = propagate!(orbp, t) + + # Compare the results. + for k = 1:length(t) + # Assemble the result vector. + st_sgp4_result = [t[k]/60 r_TEME[k]'/1000 v_TEME[k]'/1000] + + # Compare the values. + @test st_sgp4_result[1] == SGP4_results[k,1] + @test st_sgp4_result[2] ≈ SGP4_results[k,2] atol=1e-8 + @test st_sgp4_result[3] ≈ SGP4_results[k,3] atol=1e-8 + @test st_sgp4_result[4] ≈ SGP4_results[k,4] atol=1e-8 + @test st_sgp4_result[5] ≈ SGP4_results[k,5] atol=1e-9 + @test st_sgp4_result[6] ≈ SGP4_results[k,6] atol=1e-9 + @test st_sgp4_result[7] ≈ SGP4_results[k,7] atol=1e-9 + end end println("")