30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
|
# File 'lib/numru/gphys/gphys_nusdas_io.rb', line 30
def open(nusdas, varname)
if String === nusdas
nusdas = NuSDaS.new(nusdas)
end
unless NuSDaS === nusdas
raise ArgumentError, "1st arg must be a NuSDaS of path name"
end
unless (meta = nusdas.instance_variable_get(:@meta)) && meta[:nbasetime] && meta[:nbasetime] > 0
raise "NuSDaS directory is empty or nusdas_def is inconsistent with data: #{nusdas.path}"
end
var = nusdas.var(varname)
var.nil? && raise("#{varname} is not found")
data = VArrayNuSDaS.new(var)
axes = Array.new
dns = var.dim_names
lonlat2d = dns.include?("lon") && var.dim("lon").val(:full).rank > 1
dns.each{|dn|
dim = var.dim(dn)
attr = Hash.new
dim.att_names.each{|an| attr[an] = dim.att(an)}
dn << "_ref" if lonlat2d && /\A(lon|lat)\z/ =~ dn
axpos = VArray.new(dim.val(:reference), attr, dn)
axis = Axis.new(false, true)
axis.set_pos(axpos)
axes.push axis
}
grid = Grid.new(*axes)
gphys = GPhys.new(grid, data)
if lonlat2d
puts "gphys_nusdas_io.rb: assoc_coords are introduced" if $DEBUG
x = gphys.axis("lon_ref")
y = gphys.axis("lat_ref")
xy = Grid.new(x,y)
vlon = VArray.new(var.dim("lon").val(:full), x.pos.attr_copy, "lon")
vlat = VArray.new(var.dim("lat").val(:full), y.pos.attr_copy, "lat")
glon = GPhys.new(xy, vlon)
glat = GPhys.new(xy, vlat)
tmp_assoc_coords = [glon,glat]
end
if var.dim_names.include?("z") && var.att("vertical_grid") =~ /[Zz]\*/ puts "gphys_nusdas_io.rb: assoc_coords are introduced" if $DEBUG
if lonlat2d
x = gphys.axis("lon_ref")
y = gphys.axis("lat_ref")
else
x = gphys.axis("lon")
y = gphys.axis("lat")
end
z = gphys.axis("z")
xyz = Grid.new(x,y,z)
zrp, zrw, vctrans_p, vctrans_w, dvtrans_p, dvtrans_w = var.dim("z").val(:full)
zs = nusdas.var("ZSsrf") case zs.rank
when 2
na_zs = zs.get
when 3
na_zs = zs[true,true,0]
else
raise "BUG: rank of ZSsrf > 3"
end
na_zs = na_zs.newdim!(-1)
z3d = na_zs + (1 + na_zs * dvtrans_p.newdim(0,0)) * zrp.newdim(0,0)
vz3d = VArray.new(z3d, z.pos.attr_copy, "z3d")
gz3d = GPhys.new(xyz, vz3d)
tmp_assoc_coords ||= []
tmp_assoc_coords.push(gz3d)
gphys.set_att("zrp",zrp)
gphys.set_att("zrw",zrw)
gphys.set_att("vctrans_p",vctrans_p)
gphys.set_att("vctrans_w",vctrans_w)
gphys.set_att("dvtrans_p",dvtrans_p)
gphys.set_att("dvtrans_w",dvtrans_w)
end
gphys.set_assoc_coords(tmp_assoc_coords) if tmp_assoc_coords
return gphys
end
|