Granular.jl

commit 8bdb73881035e3ed26fcacbdfb4e4e134bc0eace
parent 71b3e39a5ed3f47997c736843d83c77b3d9e1a91
Author: Anders Damsgaard <andersd@riseup.net>
Date:   Thu, 22 Jun 2017 10:10:56 -0400

convert NetCDF data after read

Diffstat:
M
src/ocean.jl
 | 
33
++++++++++++++++-----------------
M
test/netcdf.jl
 | 
24
++++++++++++------------

2 files changed, 28 insertions(+), 29 deletions(-)
diff --git a/src/ocean.jl b/src/ocean.jl
@@ -89,16 +89,17 @@ function readOceanStateNetCDF(filename::String)
         error("$(filename) could not be opened")
     end
 
-    u_staggered::Array{float, 4} = NetCDF.ncread(filename, "u")
-    v_staggered::Array{float, 4} = NetCDF.ncread(filename, "v")
+    u_staggered = convert(Array{Float64, 4}, NetCDF.ncread(filename, "u"))
+    v_staggered = convert(Array{Float64, 4}, NetCDF.ncread(filename, "v"))
     u, v = interpolateOceanVelocitiesToCorners(u_staggered, v_staggered)
 
-    time = NetCDF.ncread(filename, "time")*24.*60.*60.
-    h = NetCDF.ncread(filename, "h")
-    e = NetCDF.ncread(filename, "e")
+    time = convert(Array{Float64, 1},
+                   NetCDF.ncread(filename, "time")*24.*60.*60.)
+    h = convert(Array{Float64, 4}, NetCDF.ncread(filename, "h"))
+    e = convert(Array{Float64, 4}, NetCDF.ncread(filename, "e"))
 
-    zl = NetCDF.ncread(filename, "zl")
-    zi = NetCDF.ncread(filename, "zi")
+    zl = convert(Array{Float64, 1}, NetCDF.ncread(filename, "zl"))
+    zi = convert(Array{Float64, 1}, NetCDF.ncread(filename, "zi"))
 
     return time, u, v, h, e, zl, zi
 end
@@ -110,18 +111,18 @@ Read NetCDF file with ocean *supergrid* information generated by MOM6 (e.g.
 located in the simulation `INPUT/` subdirectory.
 
 # Returns
-* `xh::Array{Float, 2}`: Longitude for cell centers [deg]
-* `yh::Array{Float, 2}`: Latitude for cell centers [deg]
-* `xq::Array{Float, 2}`: Longitude for cell corners [deg]
-* `yq::Array{Float, 2}`: Latitude for cell corners [deg]
+* `xh::Array{Float64, 2}`: Longitude for cell centers [deg]
+* `yh::Array{Float64, 2}`: Latitude for cell centers [deg]
+* `xq::Array{Float64, 2}`: Longitude for cell corners [deg]
+* `yq::Array{Float64, 2}`: Latitude for cell corners [deg]
 """
 function readOceanGridNetCDF(filename::String)
 
     if !isfile(filename)
         error("$(filename) could not be opened")
     end
-    x::Array{float, 2} = NetCDF.ncread(filename, "x")
-    y::Array{float, 2} = NetCDF.ncread(filename, "y")
+    x = convert(Array{Float64, 2}, NetCDF.ncread(filename, "x"))
+    y = convert(Array{Float64, 2}, NetCDF.ncread(filename, "y"))
 
     xh = x[2:2:end, 2:2:end]
     yh = y[2:2:end, 2:2:end]
@@ -137,16 +138,14 @@ export interpolateOceanVelocitiesToCorners
 Convert gridded data from Arakawa-C type (decomposed velocities at faces) to 
 Arakawa-B type (velocities at corners) through interpolation.
 """
-function interpolateOceanVelocitiesToCorners(u_in::Array{float, 4},
-                                             v_in::Array{float, 4})
+function interpolateOceanVelocitiesToCorners(u_in::Array{Float64, 4}, 
+                                             v_in::Array{Float64, 4})
 
     if size(u_in) != size(v_in)
         error("size of u_in ($(size(u_in))) must match v_in ($(size(v_in)))")
     end
 
     nx, ny, nz, nt = size(u_in)
-    #u = Array{float}(nx+1, ny+1, nz, nt)
-    #v = Array{float}(nx+1, ny+1, nz, nt)
     u = zeros(nx+1, ny+1, nz, nt)
     v = zeros(nx+1, ny+1, nz, nt)
     for i=1:nx
diff --git a/test/netcdf.jl b/test/netcdf.jl
@@ -28,15 +28,15 @@ u2, v2, h2, e2 = SeaIce.interpolateOceanState(ocean, ocean.time[2])
 @test_throws ErrorException SeaIce.interpolateOceanState(ocean, -1.)
 u1_5, v1_5, h1_5, e1_5 = SeaIce.interpolateOceanState(ocean,
     ocean.time[1] + (ocean.time[2] - ocean.time[1])/2.)
-@test_approx_eq u1 ocean.u[:, :, :, 1]
-@test_approx_eq v1 ocean.v[:, :, :, 1]
-@test_approx_eq h1 ocean.h[:, :, :, 1]
-@test_approx_eq e1 ocean.e[:, :, :, 1]
-@test_approx_eq u2 ocean.u[:, :, :, 2]
-@test_approx_eq v2 ocean.v[:, :, :, 2]
-@test_approx_eq h2 ocean.h[:, :, :, 2]
-@test_approx_eq e2 ocean.e[:, :, :, 2]
-@test_approx_eq u1_5 (ocean.u[:, :, :, 1] + ocean.u[:, :, :, 2])/2.
-@test_approx_eq v1_5 (ocean.v[:, :, :, 1] + ocean.v[:, :, :, 2])/2.
-@test_approx_eq h1_5 (ocean.h[:, :, :, 1] + ocean.h[:, :, :, 2])/2.
-@test_approx_eq e1_5 (ocean.e[:, :, :, 1] + ocean.e[:, :, :, 2])/2.
+@test u1 ≈ ocean.u[:, :, :, 1]
+@test v1 ≈ ocean.v[:, :, :, 1]
+@test h1 ≈ ocean.h[:, :, :, 1]
+@test e1 ≈ ocean.e[:, :, :, 1]
+@test u2 ≈ ocean.u[:, :, :, 2]
+@test v2 ≈ ocean.v[:, :, :, 2]
+@test h2 ≈ ocean.h[:, :, :, 2]
+@test e2 ≈ ocean.e[:, :, :, 2]
+@test u1_5 ≈ (ocean.u[:, :, :, 1] + ocean.u[:, :, :, 2])/2.
+@test v1_5 ≈ (ocean.v[:, :, :, 1] + ocean.v[:, :, :, 2])/2.
+@test h1_5 ≈ (ocean.h[:, :, :, 1] + ocean.h[:, :, :, 2])/2.
+@test e1_5 ≈ (ocean.e[:, :, :, 1] + ocean.e[:, :, :, 2])/2.