Merge pull request #19 from giordano/change-unit

Add method to change the unit of the result

README.md

@@ -10,13 +10,13 @@ Installation
 
 To install the package:
 
-```jlcon
+```julia
 pkg> add Cosmology
 ```
 
 Then, to load into your session:
 
-```jlcon
+```julia
 julia> using Cosmology
 ```
 
@@ -64,7 +64,7 @@ which takes the following options:
   </tr>
 </table>
 
-```jlcon
+```julia
 julia> using Cosmology
 
 julia> c = cosmology()
@@ -103,7 +103,7 @@ Distances
   </tr>
 </table>
 
-```jlcon
+```julia
 julia> using Cosmology
 
 julia> c = cosmology(OmegaM=0.26)
@@ -113,6 +113,17 @@ julia> angular_diameter_dist(c, 1.2)
 1784.0089227105113 Mpc
 ```
 
+For each function returning a unitful number, you can specify a different unit
+for the result as first argument to the function:
+
+```julia
+julia> comoving_volume(c, 0.6)
+4.936334366313069e10 Mpc^3
+
+julia> comoving_volume(u"ly^3", c, 0.6)
+1.7127035381752994e30 ly^3
+```
+
 Times
 -----
 
@@ -127,12 +138,15 @@ Times
   </tr>
 </table>
 
-```jlcon
+```julia
 julia> using Cosmology
 
 julia> c = cosmology(OmegaM=0.26)
 FlatLCDM(0.69,0.7399122024007928,0.26,8.779759920715362e-5)
 
 julia> age(c, 1.2)
 5.445600787626434 Gyr
+
+julia> lookback_time(u"yr", c, 1.2)
+8.761660748088268e9 yr
 ```

src/Cosmology.jl

@@ -2,8 +2,8 @@ __precompile__()
 
 module Cosmology
 
-using QuadGK
-using Unitful: km, s, ustrip
+using QuadGK, Unitful
+import Unitful: km, s
 using UnitfulAstro: Mpc, Gyr
 
 export cosmology,
@@ -196,6 +196,13 @@ T(c::AbstractCosmology, a0, a1; kws...) = QuadGK.quadgk(x->x/a2E(c,x), a0, a1; k
 age(c::AbstractCosmology, z; kws...) = hubble_time0(c)*T(c, 0.0, scale_factor(z); kws...)
 lookback_time(c::AbstractCosmology, z; kws...) = hubble_time0(c)*T(c, scale_factor(z), 1.0; kws...)
 
+# Easily select a different unit
+for f in (:hubble_dist0, :hubble_dist, :hubble_time0, :hubble_time, :comoving_radial_dist,
+          :comoving_transverse_dist, :angular_diameter_dist, :luminosity_dist,
+          :comoving_volume, :comoving_volume_element, :age, :lookback_time)
+    @eval $f(u::Unitful.Unitlike, args...; kws...) = uconvert(u, $f(args...; kws...))
+end
+
 ###############
 # Deprecations
 #
@@ -205,9 +212,7 @@ lookback_time(c::AbstractCosmology, z; kws...) = hubble_time0(c)*T(c, scale_fact
 @deprecate hubble_time_gyr0(c::AbstractCosmology) ustrip(hubble_time0(c::AbstractCosmology))
 @deprecate hubble_time_gyr(c::AbstractCosmology, z) ustrip(hubble_time(c::AbstractCosmology, z))
 @deprecate comoving_radial_dist_mpc(c::AbstractCosmology, z; kws...) ustrip(comoving_radial_dist(c::AbstractCosmology, z; kws...))
-@deprecate comoving_transverse_dist_mpc(c::AbstractFlatCosmology, z; kws...) ustrip(comoving_transverse_dist(c::AbstractFlatCosmology, z; kws...))
-@deprecate comoving_transverse_dist_mpc(c::AbstractOpenCosmology, z; kws...) ustrip(comoving_transverse_dist(c::AbstractOpenCosmology, z; kws...))
-@deprecate comoving_transverse_dist_mpc(c::AbstractClosedCosmology, z; kws...) ustrip(comoving_transverse_dist(c::AbstractClosedCosmology, z; kws...))
+@deprecate comoving_transverse_dist_mpc(c::AbstractCosmology, z; kws...) ustrip(comoving_transverse_dist(c::AbstractCosmology, z; kws...))
 @deprecate angular_diameter_dist_mpc(c::AbstractCosmology, z; kws...) ustrip(angular_diameter_dist(c::AbstractCosmology, z; kws...))
 @deprecate luminosity_dist_mpc(c::AbstractCosmology, z; kws...) ustrip(luminosity_dist(c::AbstractCosmology, z; kws...))
 @deprecate comoving_volume_gpc3(c::AbstractCosmology, z; kws...) ustrip(comoving_volume(c::AbstractCosmology, z; kws...))

test/runtests.jl

@@ -86,3 +86,15 @@ end
     c = cosmology(h=big(0.7), OmegaM=0.3, OmegaR=0, w0=-0.9, wa=0.1)
     @test angular_diameter_dist(c,1,rtol=dist_rtol) ≈ 1612.0585u"Mpc" rtol = dist_rtol
 end
+
+@testset "Unit conversion" begin
+    c = cosmology(h=0.9, OmegaM=0.5, OmegaR=0)
+    for u in (u"m", u"pc", u"ly")
+        @test unit(luminosity_dist(u, c, 1)) == u
+        @test unit(angular_diameter_dist(u, c, 2)) == u
+    end
+    for u in (u"s", u"yr")
+        @test unit(age(u, c, 3)) == u
+        @test unit(lookback_time(u, c, 4)) == u
+    end
+end