wcs_to_celestial_frame added

Author: aquatiko

src/Reproject.jl

@@ -3,5 +3,6 @@ module Reproject
 using FITSIO, WCS
 
 include("parsers.jl")
+include("utils.jl")
 
 end # module

src/utils.jl

@@ -0,0 +1,38 @@
+const UNDEFINED    = 987654321.0e99        # from WCS.jl
+
+"""
+    wcs_to_celestial_frame(wcs::WCSTransform)
+
+Returns the reference frame of a WCSTransform.
+The reference frame supported in Julia are FK5, ICRS and Galactic.
+"""
+function wcs_to_celestial_frame(wcs::WCSTransform)
+    radesys = wcs.radesys
+    
+    if wcs.equinox != UNDEFINED
+        equinox = wcs.equinox
+    else
+        equinox = nothing
+    end
+    
+    xcoord = wcs.ctype[1][1:4]
+    ycoord = wcs.ctype[2][1:4]
+    
+    if radesys == ""
+        if xcoord == "RA--" && ycoord == "DEC-"
+            if equinox === nothing
+                radesys = "ICRS"
+            elseif equinox < 1984.0
+                radesys = "FK4"
+            else
+                radesys = "FK5"
+            end  
+        elseif xcoord == "GLON" && ycoord == "GLAT"
+            radesys = "Gal"
+        elseif xcoord == "TLON" && ycoord == "TLAT"
+            radesys = "ITRS"
+        end
+    end
+    
+    return radesys
+end

test/runtests.jl

@@ -10,4 +10,5 @@ rp = pyimport("reproject")
 
 @testset "Reproject.jl" begin
     include("parsers.jl")
+    include("utils.jl")
 end

test/utils.jl

@@ -0,0 +1,31 @@
+@testset "wcs to celestial frame" begin
+    wcs1 = WCSTransform(2;
+                       ctype = ["RA---AIR", "DEC--AIR"],
+                       )
+    wcs2 = WCSTransform(2;
+                       ctype = ["RA---AIR", "DEC--AIR"],
+                       equinox = 1888.67 
+                       )
+    wcs3 = WCSTransform(2;
+                       ctype = ["RA---AIR", "DEC--AIR"],
+                       equinox = 2000 
+                       )
+    wcs4 = WCSTransform(2;
+                       ctype = ["GLON--", "GLAT--"], 
+                       )
+    wcs5 = WCSTransform(2;
+                       ctype = ["TLON", "TLAT"],
+                      )
+    wcs6 = WCSTransform(2;
+                       ctype = ["RA---AIR", "DEC--AIR"],
+                       radesys = "UNK" 
+                      )
+    
+    @test wcs_to_celestial_frame(wcs1) == "ICRS"
+    @test wcs_to_celestial_frame(wcs2) == "FK4"
+    @test wcs_to_celestial_frame(wcs3) == "FK5"
+    @test wcs_to_celestial_frame(wcs4) == "Gal"
+    @test wcs_to_celestial_frame(wcs5) == "ITRS"
+    @test wcs_to_celestial_frame(wcs6) == "UNK"
+end
+