Class: NumRu::GPhys
- Inherits:
-
Object
- Object
- NumRu::GPhys
- Defined in:
- lib/numru/ganalysis/fitting.rb,
lib/numru/gphys/grib.rb,
lib/numru/ganalysis/eof.rb,
lib/numru/ganalysis/met.rb,
lib/numru/gphys/ep_flux.rb,
lib/numru/gphys/version.rb,
lib/numru/gphys/gphys_io.rb,
lib/numru/ganalysis/log_p.rb,
lib/numru/gphys/gphys_fft.rb,
lib/numru/gphys/derivative.rb,
lib/numru/gphys/interpolate.rb,
lib/numru/gphys/gphys_dim_op.rb,
lib/numru/ganalysis/histogram.rb,
lib/numru/gphys/gphys_grib_io.rb,
lib/numru/ganalysis/covariance.rb,
lib/numru/gphys/coordtransform.rb,
lib/numru/gphys/gphys_grads_io.rb,
lib/numru/gphys/gphys_gtool3_io.rb,
lib/numru/gphys/gphys_io_common.rb,
lib/numru/gphys/gphys_netcdf_io.rb,
lib/numru/gphys/gphys_nusdas_io.rb,
lib/numru/gphys/gphys_hdfeos5_io.rb,
ext/numru/gphys/dim_op.c,
ext/numru/gphys/interpo.c
Overview
GPhys extension with GAnalysis::Fitting
Defined Under Namespace
Modules: Derivative, EP_Flux, GrADS_IO, GribUtils, Grib_IO, Gtool3_IO, HE5_IO, IO, IO_Common, NetCDF_IO, NuSDaS_IO Classes: Grib, GribDim, GribVar
Constant Summary collapse
- VERSION =
Add alpha while under development; remove it to release
"1.5.0"
- COS_TAPER_SP_FACTOR =
1.0 / 0.875
- BC_SIMPLE =
enum in convol_filter.c
10
- BC_CYCLIC =
enum in convol_filter.c
11
- BC_TRIM =
enum in convol_filter.c
12
- @@fft_forward =
-1
- @@fft_backward =
1
- @@fft_ignore_missing =
false
- @@fft_missing_replace_val =
nil
- @@interpo_previous_cutter =
nil
- @@interpo_previous_modifier =
nil
- @@interpo_missval =
NC_FILL_DOUBLE/FLOAT ~15*2^119
9.9692099683868690e+36
- @@interpo_extrapolation =
false
- @@default_missval =
NC_FILL_DOUBLE/FLOAT ~15*2^119
9.9692099683868690e+36
Class Method Summary collapse
-
.c_cap_by_boundary ⇒ Object
cap_by_boundary : Cap (insert) a NArray with boundary values.
-
.c_cell_integ_irreg ⇒ Object
cell_integ_irreg: trapezoidal numerical integration over coordinate cells, supporting irregular grid.
-
.c_cum_integ_irreg ⇒ Object
cum_integ_irreg : similar to cell_integ_irreg but it sums up along the axis.
-
.extrapolation=(extrapo) ⇒ Object
Change the behavior of the interpolation methods to extrapolate outside the grid coverage.
- .fft_ignore_missing(ignore = true, replace_val = nil) ⇒ Object
-
.interpo_find_loc_1D ⇒ Object
to make “find loc” methods available outside GPhys as class methods.
-
.interpo_find_loc_1D_MD ⇒ Object
To apply interpo_find_loc_1D multi-dimensionally.
Instance Method Summary collapse
-
#bin_mean(dim, len, nminvalid = 1) ⇒ Object
Binning along a dimension (mean).
-
#bin_sum(dim, len, nminvalid = 1) ⇒ Object
Binning along a dimension (summation).
- #cderiv(*args) ⇒ Object
-
#coord_data_reverse(axname, pos) ⇒ Object
Reverse the main data (i.e., the dependent variable) and one of the coordinates (an independent variable) through interpolation.
- #coordtransform(coordmapping, axes_to, *dims) ⇒ Object
- #corelation(other, *dims) ⇒ Object (also: #correlation)
-
#cos_taper(*dims) ⇒ Object
Spectral factor for the cosine taper.
- #covariance(other, *dims) ⇒ Object
-
#dcl_fig_cut(dimx, dimy, ux, uy) ⇒ Object
Interpolation on the DCL window (automatic iso-interval interpolation along a poly line that can be drawn in the current viewport of the DCL window).
- #deriv2nd(*args) ⇒ Object
- #detrend(*dims) ⇒ Object
- #eof(*args) ⇒ Object
- #fft(backward = false, *dims) ⇒ Object
- #fft_deriv(dim) ⇒ Object
- #histogram(opts = Hash.new) ⇒ Object (also: #histogram1D)
-
#interpolate(*coords) ⇒ Object
Wide-purpose multi-dimensional linear interpolation.
- #least_square_fit(functions, ensemble_dims = nil, indep_dims = nil) ⇒ Object
-
#logp_coord_p2z(pdim = nil) ⇒ Object
Convert the pressure coordinate in self to log-pressure height (after duplicating self).
-
#mouse_cut(dimx, dimy, num = 2, line_type = 1, line_index = 1) ⇒ Object
Makes a subset interactively by specifying a (poly-)line on the DCL viewport.
-
#mouse_cut_repeat ⇒ Object
Interpolation onto grid points specified by the previous call of GPhys#mouse_cut.
- #phase_velocity(kdim, fdim, kconv, fconv, kf0_is_c0 = true, no_kfreorder = false) ⇒ Object
- #phase_velocity_binning(kdim, fdim, cbins, kconv = nil, fconv = nil) ⇒ Object
- #phase_velocity_binning_iso_norml(kdim, fdim, cmin, cmax, cint, kconv = nil, fconv = nil) ⇒ Object
- #phase_velocity_filter(xdim, tdim, cmin = nil, cmax = nil, xconv = nil, tconv = nil, remove_xtmean = false) ⇒ Object
- #rawspect2powerspect(*dims) ⇒ Object
-
#regrid(to) ⇒ Object
Interpolate to conform the grid to a target GPhys object.
-
#running_mean(dim, len_or_wgt = nil, bc = BC_SIMPLE, nminvalid = 1) ⇒ Object
Running mean along a dimension (with optional weight specification).
- #spect_one_sided(dim) ⇒ Object
- #spect_zero_centering(dim) ⇒ Object
- #threepoint_O2nd_deriv(*args) ⇒ Object
Class Method Details
.c_cap_by_boundary ⇒ Object
cap_by_boundary : Cap (insert) a NArray with boundary values
Restriction; data alignment is restricted so that the beginning of the out data is always valid (within the domain). To ensure it, it should be either zcrd is increasing and upper==true or zcrd is decreasing and upper==false.
RETURN VALUES
fe: f capped by the boundary values. The dimension zdim is
extended by 1; i.e., f[:,nz,:] --> fe[:,nz+1,:], where ":" respresent
arbitrary number of dimensions. The elements of fe are equal to
those of f where inside the domain (simple copies), and they are equal
to the elements of fs at the bondary (simple copies if fs is given;
if not, guessed by interpolation or naive extension).
ze: grid points of fe along zdim. It is a mixture of zcrd and zs;
it is zcrd inside the domain (where f is copied), and it is zs
at the boundary (where fs is copied).
Same shape as fe.
nze: The number of valid data along zdim of fe. Shaped as ze,
according to the notation above. For example, when fe is 4D and
zdim==2, fe[i,j,k,l] is valid for k = 0,1,...,nze[i,j,l]-1,
where the boundary is at nze[i,j,l]-1. Thus, nze is always
smaller than or equal to the length of zdim of fe (which is nz+1)
1023 1024 1025 |
# File 'ext/numru/gphys/dim_op.c', line 1023 static VALUE cap_by_boundary(obj, f, zdim, zcrd, upper, zb, fb) VALUE obj; |
.c_cell_integ_irreg ⇒ Object
cell_integ_irreg: trapezoidal numerical integration over coordinate cells, supporting irregular grid
Description
Suppose a multi-dimensional NArray f, where colon represents any number of dimensions, and k is the “z” dimension along which integration is made. We write its real space representation as f(z; x), where x symbolically represents all of the independent variables other than z, and for simplicity, we further write it as f(z).
z is sampled at z_k, k=0,1,…,nzbound-1. This method allows z_k to be defined for each z column, so it requires a multi-D NArray argument z (having the same shape as f). Optionally, nzbound can also vary as nzbound. If, instead, nil is given to nzbound, the entire z grid is used; nzbound is set to f.shape(zdim).
We define the integration of f as
{ \int_za^zb f(z) dz, when za<=zb,
I(za,zb) = {
{ -\int_za^zb f(z) dz, otherwise.
In other words, our integration is always made from the smaller end to the greater end.
In the normal use case (when w is given nil), we define the cell integration as,
I(-\infty, zc_0), I(zc_0, zc_1), I(zc_1, zc_2),...,
The cell boundaries zc_m (m=0,1,..) are specified by the 1D NArray argument “ccell”; ccell must be aligned in the increasing order.
This method allows coordinate transformation by specifying another coordinate variable w (having the same shape as f). In this case, the ccell argument specifies a coordinate with respect to w: wc_m (m=0,1,…; wc_m must be in the increasing order). The integration is still taken with respect to z, so the cell integration is expressed as
I(-\infty, z(wc_0)), I(z(wc_0), z(wc_1)), I(z(wc_1), z(wc_2)),...,
The grid values z and w do not have to be monotonic; the numerical integration properly treats the contribution from multiple ranges along k. Mathematically, the coordinate-transferred integration over the w bin (-infty, wc] is expressed as
\int_-\infty^+\infty H(wc-w(z)) f(z) dz,
where H is the Heaviside function. The normal use case (without w) is simply when w is z itself, which is exploited in implementation.
790 791 792 |
# File 'ext/numru/gphys/dim_op.c', line 790 static VALUE cell_integ_irreg(obj, f, z, zdim, nzbound, ccell, w) VALUE obj; |
.c_cum_integ_irreg ⇒ Object
cum_integ_irreg : similar to cell_integ_irreg but it sums up along the axis. – This method acutually uses cell_integ_irreg and make sumation.
979 980 981 |
# File 'ext/numru/gphys/dim_op.c', line 979 static VALUE cum_integ_irreg(obj, f, z, zdim, nzbound, ccell, w) VALUE obj; |
.extrapolation=(extrapo) ⇒ Object
Change the behavior of the interpolation methods to extrapolate outside the grid coverage.
ARGUMENTS
-
extrapo : true or false — the default behaviour is false (not to extrapolate), so use this method if you want to set it to true.
28 29 30 |
# File 'lib/numru/gphys/interpolate.rb', line 28 def self.extrapolation=(extrapo) @@interpo_extrapolation = extrapo end |
.fft_ignore_missing(ignore = true, replace_val = nil) ⇒ Object
331 332 333 334 |
# File 'lib/numru/gphys/gphys_fft.rb', line 331 def self.fft_ignore_missing( ignore=true, replace_val=nil ) @@fft_ignore_missing = ignore @@fft_missing_replace_val = replace_val end |
.interpo_find_loc_1D ⇒ Object
to make “find loc” methods available outside GPhys as class methods
311 312 313 |
# File 'ext/numru/gphys/interpo.c', line 311 static VALUE interpo_find_loc_1D(obj, X, x, missval, extrapo) VALUE obj; |
.interpo_find_loc_1D_MD ⇒ Object
To apply interpo_find_loc_1D multi-dimensionally
372 373 374 |
# File 'ext/numru/gphys/interpo.c', line 372 static VALUE interpo_find_loc_1D_MD(obj, X, x, dimc, missval, extrapo) VALUE obj; |
Instance Method Details
#bin_mean(dim, len, nminvalid = 1) ⇒ Object
Binning along a dimension (mean)
The values are averaged every “len” grids; unlike running_mean the number of grids is reduced to 1/len. Currently, the only available boundary condition is BC_TRIM.
ARGUMENTS
-
dim (Integer or String) : the dimension
-
len (Integer): length of the bin
-
nminvalid (Integer; optional; defualt=1): Effective only for data with missing. Minimum number of grid points needed for averaging (1 to len).
RETURN VALUE
-
a GPhys
113 114 115 116 |
# File 'lib/numru/gphys/gphys_dim_op.rb', line 113 def bin_mean(dim, len, nminvalid=1) dim = dim_index(dim) # to handle String or negative specification GPhys.new( grid.binning(dim, len), data.bin_mean(dim, len, nminvalid) ) end |
#bin_sum(dim, len, nminvalid = 1) ⇒ Object
Binning along a dimension (summation)
Similar to bin_mean, but the values are simply summed without averaging
ARGUMENTS
-
dim (Integer or String) : the dimension
-
len (Integer): length of the bin
-
nminvalid (Integer; optional; defualt=1): Effective only for data with missing. Minimum number of grid points needed for averaging (1 to len).
RETURN VALUE
-
a GPhys
131 132 133 134 |
# File 'lib/numru/gphys/gphys_dim_op.rb', line 131 def bin_sum(dim, len, nminvalid=1) dim = dim_index(dim) # to handle String or negative specification GPhys.new( grid.binning(dim, len), data.bin_sum(dim, len, nminvalid) ) end |
#cderiv(*args) ⇒ Object
91 92 93 |
# File 'lib/numru/gphys/derivative.rb', line 91 def cderiv(*args) Derivative::cderiv(self,*args) end |
#coord_data_reverse(axname, pos) ⇒ Object
Reverse the main data (i.e., the dependent variable) and one of the coordinates (an independent variable) through interpolation.
Returns a GPhys in which the main data is the specfied coordinate (argument: axname) sampled at specified locations (argument: pos) in terms of the main data of self. The main data of self is expected to be quai-monotonic with respect to the specfied coordinate.
ARGUMENTS
-
axname [String] : one of the names of the axes (i.e. main coordinates. Auxiliary coordinates are not supported as the target.)
-
pos [NArray] : grid locations. For example, if the current data is potential temperature theta, pos consists of the theta levels to make sampling.
RETURN VALUE
-
a GPhys
231 232 233 234 235 236 237 238 239 240 |
# File 'lib/numru/gphys/interpolate.rb', line 231 def coord_data_reverse(axname,pos) gp = self.axis(axname).to_gphys gp = self.shape_coerce_full(gp)[0] # conform the shape to that of self gp = GPhys.new( gp.grid.copy, gp.data ) # copy grid to avoid side effect # on the grid of self gp.set_assoc_coords([self]) pos = NArray[*pos].to_type(NArray::FLOAT) if pos.is_a?(Array) newcrd = VArray.new(pos,self.data,self.name) # succeeds the attributes gp.interpolate(axname=>newcrd) end |
#coordtransform(coordmapping, axes_to, *dims) ⇒ Object
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 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 |
# File 'lib/numru/gphys/coordtransform.rb', line 41 def coordtransform( coordmapping, axes_to, *dims ) rankmp = coordmapping.rank #< check arguments > if axes_to.length != rankmp raise ArgumentError, "length of axes_to must be equal to the rank of coordmapping" end if self.rank == rankmp dims = (0...rankmp).collect{|i| i} elsif self.rank < rankmp raise ArgumentError,"rank of coordmapping is greater than self.rank" elsif dims.length != rankmp raise ArguemntError, "# of dimensions speficied is not equal to the rank of coordmapping" elsif dims != dims.sort raise ArguementErroor,"dims must be in the increasing order" end #< get grid points > vt = coordmapping.map_grid( *dims.collect{|d| axes_to[d].pos.val} ) x = dims.collect{|d| self.grid.axis(d).pos.val} #< prepare the output object > axes = (0...self.rank).collect{|i| grid.axis(i)} dims.each_with_index{|d,j| axes[d]=axes_to[j]} grid_to = Grid.new( *axes ) vnew = VArray.new( NArray.new( self.data.ntype, *grid_to.shape ), self.data, self.name ) #< do interpolation (so far only 2D is supported) > case dims.length when 2 if !HAVE_NUMRU_SSL2 p "interpolation without SSL2" # raise "Sorry, so far I need SSL2 (ruby-ssl2)" self.( *dims ){ |fxy,idx| wgts = Array.new idxs = Array.new for d in 0..dims.length-1 wgt = vt[d].dup.fill!(-1.0) idx0 = vt[d].dup.to_i.fill!(-1) idx1 = idx0.dup.fill!(x[d].length) xsort = x[d].sort xsortindex = x[d].sort_index for i in 0..x[d].length-1 idx0[ xsort[i] <= vt[d] ] = xsortindex[i] idx1[ xsort[-1-i] >= vt[d] ] = xsortindex[-1-i] end # where idx0=idx1 wgt[ idx0.eq(idx1) ] = 1.0 # where vt[d] < x[d].min wgt[ idx0 <= -1 ] = 1.0 idx0[ idx0 <= -1 ] = 0 # where vt[d] > x[d].max wgt[ idx1 >= x[d].length ] = 0.0 idx1[ idx1 >= x[d].length ] = x[d].length-1 # normal points mask = wgt.eq(-1.0) wgt[mask] = (vt[d][mask]-x[d][idx0[mask]])/(x[d][idx1[mask]]-x[d][idx0[mask]]) wgts.push(wgt) idxs[d*2] = idx0 idxs[d*2+1] = idx1 end case dims.length # when 1 # f = fxy.data.val[idxs[0]]*(1-wgts[0]) + # fxy.data.val[idxs[1]]*wgts[0] # f = f.to_na if( f.class.to_s == "NArrayMiss" ) when 2 lx = fxy.shape[0] f = ( fxy.data.val[idxs[0]+idxs[2]*lx]*(1-wgts[0]) + fxy.data.val[idxs[1]+idxs[2]*lx]*wgts[0] ) * (1-wgts[1]) + ( fxy.data.val[idxs[0]+idxs[3]*lx]*(1-wgts[0]) + fxy.data.val[idxs[1]+idxs[3]*lx]*wgts[0] ) * wgts[1] f = f.to_na if( f.class.to_s == "NArrayMiss" ) else raise "Sorry, #{v.length}D interpolation is yet to be supported" end if(idx==false) vnew[] = f else vnew[*idx] = f end } else ix=iy=0 m=3 self.( *dims ){ |fxy,idx| c,xt = SSL2.bicd3(x[0],x[1],fxy.val,m) begin ix,iy,f = SSL2.bifd3(x[0],x[1],m,c,xt,0,vt[0],ix,0,vt[1],iy) rescue $stderr.print "Interpolation into", vt[0].inspect, vt[1].inspect raise $! end vnew[*idx] = f } end else raise "Sorry, #{v.length}D interpolation is yet to be supported" end #< finish > GPhys.new( grid_to, vnew ) end |
#corelation(other, *dims) ⇒ Object Also known as: correlation
93 94 95 |
# File 'lib/numru/ganalysis/covariance.rb', line 93 def corelation(other, *dims) GAnalysis.corelation(self, other, *dims) end |
#cos_taper(*dims) ⇒ Object
Spectral factor for the cosine taper. Specta should be multiplied by this.
339 340 341 342 343 344 345 346 347 348 349 350 351 352 353 354 355 356 357 358 359 360 361 362 |
# File 'lib/numru/gphys/gphys_fft.rb', line 339 def cos_taper(*dims) if dims.length < 1 raise ArgumentError,'You have to specify one or more dimensions' end dims.sort!.uniq! val = self.data.val dims.each{|dim| dim = dim_index(dim) if dim.is_a?(String) dim += rank if dim < 0 raise ArgumentError,"dim #{dim} does not exist" if dim<0 || dim>rank nx = shape[dim] wgt = NArray.float(nx).fill!(1) x = 10.0 / nx * (NArray.float(nx).indgen!+0.5) wskl = x.lt(1).where wskr = x.gt(9).where wgt[wskl] = 0.5*( 1.0 - NMath::cos(Math::PI*x[wskl]) ) wgt[wskr] = 0.5*( 1.0 - NMath::cos(Math::PI*x[wskr]) ) wgt.reshape!( *([1]*dim + [nx] + [1]*(rank-dim-1)) ) val = val*wgt } to_ret = self.copy to_ret.data.val = val to_ret end |
#covariance(other, *dims) ⇒ Object
89 90 91 |
# File 'lib/numru/ganalysis/covariance.rb', line 89 def covariance(other, *dims) GAnalysis.covariance(self, other, *dims) end |
#dcl_fig_cut(dimx, dimy, ux, uy) ⇒ Object
Interpolation on the DCL window (automatic iso-interval interpolation along a poly line that can be drawn in the current viewport of the DCL window). Used in mouse_cut.
ARGUMENTS
-
dimx [Integer or String] : specifies the dimension corresponding to the UX coordinate. (Here, the UX coordinate is the X coordinate of the DCL’s USER coordinate. For exapmle, longitude if map projection.)
-
dimy [Integer or String] : specifies the dimension corresponding to the UY coordinate. (Here, the UY coordinate is the Y coordinate of the DCL’s USER coordinate. For exapmle, latitude if map projection.)
-
ux [Array] : x values in terms of the UX coordinate
-
uy [Array] : y values in terms of the UY coordinate Lengths of ux and uy must be the same and greter or equal to 2.
89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 |
# File 'lib/numru/gphys/interpolate.rb', line 89 def dcl_fig_cut(dimx,dimy,ux,uy) len = ux.length raise("ux and uy must be arrays with the (same) length >= 2") if len<=1 raise("ux's len (#{len}) != uy's len (#{uy.length})") if uy.length != len vx=Array.new; vy=Array.new for i in 0...len vx[i],vy[i] = NumRu::DCL.stftrf(ux[i],uy[i]) end kx = Array.new ky = Array.new cut = [true]*rank for i in 0...len cut[dimx] = ux[i] cut[dimy] = uy[i] dummy, sl = grid.cut(*cut) kx[i] = sl[dimx] ky[i] = sl[dimy] end ndiv = Array.new ndsum = [0] for i in 0...len-1 ndiv[i] = Math.sqrt( (kx[i+1]-kx[i])**2 + (ky[i+1]-ky[i])**2).to_i ndiv[i] += 1 if i==len-2 ndsum.push ndsum[-1] + ndiv[i] # 0, ndiv[0], ndiv[0]+ndiv[1], ... end ndtot = ndsum[-1] vxdiv = NArray.float(ndtot) vydiv = NArray.float(ndtot) for i in 0...len-1 if i!=len-2 a = NArray.float(ndiv[i]).indgen / ndiv[i] else a = NArray.float(ndiv[i]).indgen / (ndiv[i]-1) end vxdiv[ndsum[i]...ndsum[i+1]] = (1.0-a)*vx[i] + a*vx[i+1] vydiv[ndsum[i]...ndsum[i+1]] = (1.0-a)*vy[i] + a*vy[i+1] end uxdiv = NArray.float(ndtot) uydiv = NArray.float(ndtot) for i in 0...ndtot uxdiv[i], uydiv[i] = DCL.stitrf(vxdiv[i], vydiv[i]) end cx = coord(dimx) xcrd = VArray.new(uxdiv, cx, cx.name) cy = coord(dimy) ycrd = VArray.new(uydiv, cy, cy.name) if (vxdiv[-1]-vxdiv[0]).abs > (vydiv[-1]-vydiv[0]).abs cutter = [xcrd,ycrd] # x will be the main coord var if not map proj crd = xcrd else cutter = [ycrd,xcrd] # x will be the main coord var if not map proj crd = ycrd end axnm = crd.name itr = DCL.sgqtrn if itr>=10 and itr<=40 newcrd = __sp_dist(xcrd,ycrd) modifier = Proc.new{|gp| newax = Axis.new.set_pos(newcrd) gp.grid.set_axis(axnm,newax) g = Grid.new( newax ) gxcrd = GPhys.new(g,xcrd) gycrd = GPhys.new(g,ycrd) gp.set_assoc_coords([gxcrd, gycrd]) gp } else modifier = nil end @@interpo_previous_cutter = cutter @@interpo_previous_modifier = modifier # < do the job > gpnew = interpolate(cutter) gpnew = modifier[gpnew] if modifier gpnew end |
#deriv2nd(*args) ⇒ Object
97 98 99 |
# File 'lib/numru/gphys/derivative.rb', line 97 def deriv2nd(*args) Derivative::deriv2nd(self,*args) end |
#detrend(*dims) ⇒ Object
364 365 366 367 368 369 370 371 372 373 374 375 376 377 378 379 380 381 382 383 384 385 386 387 388 389 390 391 392 393 394 395 396 397 398 399 400 401 402 403 404 405 406 407 408 409 410 411 412 413 414 |
# File 'lib/numru/gphys/gphys_fft.rb', line 364 def detrend(*dims) if dims.length < 1 raise ArgumentError,'You have to specify one or more dimensions' end dims.sort!.uniq! val = self.data.val dims.each{|dim| dim = dim_index(dim) if dim.is_a?(String) dim += rank if dim < 0 raise ArgumentError,"dim #{dim} does not exist" if dim<0 || dim>rank if val.is_a?(NArray) x = self.coord(dim).val x.reshape!( *([1]*dim + [x.length] + [1]*(rank-dim-1)) ) vmean = val.mean(dim) vxmean = (val*x).mean(dim) xmean = x.mean(dim) x2mean = (x*x).mean(dim) denom = x2mean-xmean**2 if denom != 0 a = (vxmean - vmean*xmean)/denom b = (vmean*x2mean - vxmean*xmean)/denom else a = 0 b = vmean end elsif val.is_a?(NArrayMiss) x = self.coord(dim).val x.reshape!( *([1]*dim + [x.length] + [1]*(rank-dim-1)) ) x = NArrayMiss.to_nam( NArray.new(x.typecode, *val.shape) + x, val.get_mask ) vmean = val.mean(dim) vxmean = (val*x).mean(dim) xmean = x.mean(dim) x2mean = (x*x).mean(dim) denom = x2mean-xmean**2 meq0 = denom.eq(0).to_na(0) # ==0 and not masked mne0 = denom.ne(0).to_na(0) # !=0 and not masked denom.set_mask(mne0) # only nonzero part will be used to divide: a = (vxmean - vmean*xmean)/denom b = (vmean*x2mean - vxmean*xmean)/denom a[meq0] = 0 b[meq0] = vmean[meq0] end a.newdim!(dim) if !a.is_a?(Numeric) b.newdim!(dim) if !b.is_a?(Numeric) val = val - a*x-b } to_ret = self.copy to_ret.data.val = val to_ret end |
#eof(*args) ⇒ Object
229 230 231 |
# File 'lib/numru/ganalysis/eof.rb', line 229 def eof(*args) GAnalysis.eof(self, *args) end |
#fft(backward = false, *dims) ⇒ Object
416 417 418 419 420 421 422 423 424 425 426 427 428 429 430 431 432 433 434 435 436 437 438 439 440 441 442 443 444 445 446 447 448 449 450 451 452 453 454 455 456 457 458 459 460 461 462 463 464 465 466 467 468 469 470 471 472 473 474 475 476 477 478 479 480 481 482 483 484 485 486 487 488 489 490 491 492 493 494 495 496 497 498 499 500 501 502 503 |
# File 'lib/numru/gphys/gphys_fft.rb', line 416 def fft(backward=false, *dims) fftw3 = false if defined?(FFTW3) fftw3 = true elsif !defined?(FFTW) raise "Both FFTW3 and FFTW are not installed." end if backward==true dir = @@fft_backward elsif !backward dir = @@fft_forward else raise ArgumentError,"1st arg must be true or false (or, equivalenty, nil)" end # <FFT> gfc = self.copy # make a deep clone if fftw3 val = gfc.data.val if @@fft_ignore_missing and val.is_a?(NArrayMiss) if @@fft_missing_replace_val val = val.to_na(@@fft_missing_replace_val) else val = val.to_na end elsif val.is_a?(NArrayMiss) && val.count_invalid == 0 val = val.to_na end fcoef = FFTW3.fft( val, dir, *dims ) else # --> always FFT for all dimensions if dims.length == 0 raise ArgumentError, "dimension specification is available only if FFTW3 is installed" end val = gfc.data.val if @@fft_ignore_missing and val.is_a?(NArrayMiss) if @@fft_missing_replace_val val = val.to_na(@@fft_missing_replace_val) else val = val.to_na end elsif val.is_a?(NArrayMiss) && val.count_invalid == 0 val = val.to_na end fcoef = FFTW.fftw( val, dir ) end if dir == @@fft_forward if dims.length == 0 fcoef = fcoef / fcoef.length.to_f # normalized if forward FT else sh = fcoef.shape len = 1 dims.each{|d| raise ArgumentError, "dimension out of range" if sh[d] == nil len *= sh[d] } fcoef = fcoef / len end end gfc.data.replace_val( fcoef ) # <coordinate variables> for i in 0...gfc.rank if dims.length == 0 || dims.include?(i) || dims.include?(i+rank) __predefined_coord_units_conversion(gfc.coord(i)) cv = gfc.coord(i).val n = cv.length clen = (cv.max - cv.min) * n / (n-1) wn = (2*Math::PI/clen) * NArray.new(cv.typecode,cv.length).indgen!.to_f if (!backward) gfc.coord(i).set_att('origin_in_real_space',cv[0..0]) else if ( org = gfc.coord(i).get_att('origin_in_real_space') ) wn += org[0] ###gfc.coord(i).del_att('origin_in_real_space') end end gfc.coord(i).replace_val(wn) gfc.coord(i).units = gfc.coord(i).units**(-1) __coord_name_conversion(gfc.coord(i), backward) end end # <fini> gfc end |
#fft_deriv(dim) ⇒ Object
505 506 507 508 509 510 511 512 513 514 515 516 517 518 519 520 521 522 523 524 525 526 527 |
# File 'lib/numru/gphys/gphys_fft.rb', line 505 def fft_deriv(dim) tp = self.data.typecode fc = self.fft(false,dim) wn = fc.coord(dim) k = wn.val.to_type(NArray::Complex) n = k.length n2a = (n-1)/2 n2b = [n/2 + 1, n-1].min # min to avoid error if n=2 (though meaningless) kmx = k[-1]+k[1] ik = NArray.complex(n) ik[0..n2a] = k[0..n2a]*Complex::I ik[n2b..-1] = (k[n2b..-1]-kmx) * Complex::I dim.times{ik.newdim!(0)} (self.rank-dim-1).times{ik.newdim!(-1)} fc.replace_val(fc.val*ik) deriv = fc.fft(true,dim) deriv.units = deriv.units * wn.units if tp >= NArray::SCOMPLEX deriv else deriv.real end end |
#histogram(opts = Hash.new) ⇒ Object Also known as: histogram1D
155 156 157 |
# File 'lib/numru/ganalysis/histogram.rb', line 155 def histogram(opts=Hash.new) GAnalysis.histogram(self, opts) end |
#interpolate(*coords) ⇒ Object
Wide-purpose multi-dimensional linear interpolation
This method supports interpolation regarding combinations of 1D and 2D coordinate variables. For instance, suppose self is 4D with coordinates named [“x”, “y”, “z”, “t”] and associated coordinates “sigma” (“sigma” is 1D and its axis is “z”), “p”, “q” (“p” and “q” are 2D having the coordinates “x” and “y”). You can make interpolation by specifying 1D VArrays whose names are among “x”, “y”, “z”, “t”, “sigma”, “p”, “q”. You can also use a Hash like => 1.0 to specify a single point along the “x” coordinate.
If the units of the target coordinate and the current coordinate are different, a converstion was made so that slicing is made correctly, as long as the two units are comvertible; if the units are not convertible, it is just warned.
If you specify only “x”, “y”, and “t” coordinates for interpolation, the remaining coordinates “z” is simply retained. So the result will be 4 dimensional with coordinates named [“x”, “y”, “z”, “t”], but the lengths of “x”, “y”, and “t” dimensions are changed according to the specification. Note that the result could be 3-or-smaller dimensional – see below.
Suppose you have two 1D VArrays, xnew and ynew, having names “x” and “y”, respectively, and the lengths of xnew and the ynew are the same. Then, you can give an array of the two, [xnew, ynew], for coord0 as
gp_int = gp_org.interpolate( [xnew, ynew] )
(Here, gp_org represents a GPhys object, and the return value pointed by gp_int is also a GPhys.) In this case, the 1st dimension of the result (gp_int) will be sampled at the points [xnew,ynew], [xnew,ynew], [xnew,ynew], …, while the 2nd and the third dimensions are “z” and “t” (no interpolation). This way, the rank of the result will be reduced from that of self.
If you instead give xnew to coord0 and ynew to coord1 as
gp_int = gp_org.interpolate( xnew, ynew )
The result will be 4-dimensional with the first coordinate sampled at xnew, xnew, xnew,… and the second coordinate sampled at ynew, ynew, ynew,…
You can also cut regarding 2D coordinate variable as
gp_int = gp_org.interpolate( pnew, qnew )
gp_int = gp_org.interpolate( xnew, qnew )
gp_int = gp_org.interpolate( [pnew, qnew] )
gp_int = gp_org.interpolate( [xnew, qnew] )
In any case, the desitination VArrays such as xnew ynew pnew qnew must be one-dimensional.
Note that
gp_int = gp_org.interpolate( qnew )
fails (exception raised), since it is ambiguous. If you tempted to do so, perhaps what you want is covered by the following special form:
As a special form, you can specify a particular dimension like this:
gp_int = gp_org.interpolate( "x"=>pnew )
Here, interpolation along “x” is made, while other axes are retained. This is useful if pnew corresponds to a multi-D coordinate variable where there are two or more corresponding axes (otherwise, this special form is not needed.)
See the test part at the end of this file for more examples.
LIMITATION
Currently associated coordinates expressed by 3D or greater dimensional arrays are not supported.
Computational efficiency of pure two-dimensional coordinate support should be improved by letting C extensions cover deeper and improving the search algorithm for grid (which is usually ordered quasi-regularly).
COVERAGE
Extrapolation is covered for 1D coordinates, but only interpolation is covered for 2D coordinates (which is limited by gt2dlib in DCL – exception will be raised if you specify a grid point outside the original 2D grid points.).
MATHEMATICAL SPECIFICATION
The multi-dimensional linear interpolation is done by supposing a (hyper-) “rectangular” grid, where each dimension is independently sampled one-dimensionally. In case of interpolation along two dimensional coordinates such as “p” and “q” in the example above, a mapping from a rectangular grid is assumed, and the corresponding points in the rectangular grid is solved inversely (currently by using gt2dlib in DCL).
For 1D and 2D cases, linear interpolations may be expressed as
1D: zi = (1-a)*z0 + a*z1
2D: zi = (1-a)*(1-b)*z00 + a*(1-b)*z10 + (1-a)*b*z01 + a*b*z11
This method is extended to arbitrary number of dimensions. Thus, if the number of dimensions to interpolate is S, then 2**S grid points are used for each interpolation (8 points for 3D, 16 points for 4D,…). Thus, the linearity of this interpolation is only along each dimension, not over the whole dimensionality.
USAGE
interpolate(coord0, coord1, ...)
ARGUMENTS
-
coord0, coord1,… [ 1D VArray, or Array of 1D VArray, or a 1-element Hash as => slice_loc_value(Numeric) ] : locations to which interpolation is made. Names of all the VArray’s in the arguments must exist among the names of the coordinates of self (including associated coordinates), since the dimension finding is made in terms of coordinate names. If an argument is an Array of VArray’s, the first VArray will become the main coordinate variable, and the rest will be associated coordinates.
- SPECIAL CASE
-
You can specfify a one-element Hash as the only argument such as
gphys.interpolate("x"=>varray)
where varray is a coordinate onto which interpolation is made. This is espcially useful if varray is multi-D. If varray’s name “p” (name of a 2D coordnate var), for example, you can interpolate only regarding “x” by retaining other axes. If varray is 1-diemnsional, the same thing can be done simply by
gphys.interpolate(varray)
since the corresponding 1D coordinate is found aotomatically.
RETURN VALUE
-
a GPhys
389 390 391 392 393 394 395 396 397 398 399 400 401 402 403 404 405 406 407 408 409 410 411 412 413 414 415 416 417 418 419 420 421 422 423 424 425 |
# File 'lib/numru/gphys/interpolate.rb', line 389 def interpolate(*coords) coords, org_coords, org_dims, newgrid = _interpo_match_coords(coords) crdmap = _interpo_reorder_2crdmap(coords, org_coords, org_dims) idxmap = _interpo_find_position(crdmap) z = val if z.is_a?(NArrayMiss) missval = ( (a=get_att('_FillValue')) ? a[0] : nil ) || ( (a=get_att('missing_value')) ? a[0] : nil ) || @@interpo_missval z = z.to_na(missval) input_nomiss = false else input_nomiss = true if @@interpo_extrapolation missval = nil else missval = @@interpo_missval end end na = c_interpo_do(newgrid.shape, idxmap, z, missval, @@interpo_extrapolation) # [C-extension] if !input_nomiss || !@@interpo_extrapolation mask = na.ne(missval) if !input_nomiss || mask.min == 0 na = NArrayMiss.to_nam_no_dup(na,mask) end end va = VArray.new(na, data, name) ret = GPhys.new(newgrid, va) ret.grid.set_lost_axes(self.lost_axes) ret end |
#least_square_fit(functions, ensemble_dims = nil, indep_dims = nil) ⇒ Object
460 461 462 463 464 465 466 467 468 469 470 471 472 473 474 475 476 477 478 479 480 481 482 483 484 485 486 487 488 489 490 491 492 493 494 495 496 497 498 499 500 |
# File 'lib/numru/ganalysis/fitting.rb', line 460 def least_square_fit(functions, ensemble_dims=nil, indep_dims=nil) #< preparation > no_fitting_dims = Array.new if ensemble_dims ensemble_dims = ensemble_dims.collect{|d| @grid.dim_index(d)} no_fitting_dims += ensemble_dims end if indep_dims indep_dims = indep_dims.collect{|d| @grid.dim_index(d)} no_fitting_dims += indep_dims end fitting_dims = (0...rank).collect{|i| i} - no_fitting_dims grid_locs = fitting_dims.collect{|d| coord(d).val} data = self.val #< fitting > c, bf, diff = GAnalysis::Fitting.least_square_fit(data, grid_locs, functions, ensemble_dims, indep_dims) #< make a GPhys of the best fit > if !ensemble_dims grid = self.grid else axes = Array.new (0...rank).each{|d| axes.push(self.axis(d)) unless ensemble_dims.include?(d) } grid = Grid.new(*axes) shape = bf.shape ensemble_dims.sort.reverse_each{|d| shape.delete_at(d)} bf = bf.reshape(*shape) end va = VArray.new(bf, self.data, self.name) bf = GPhys.new(grid, va) [c, bf, diff] end |
#logp_coord_p2z(pdim = nil) ⇒ Object
Convert the pressure coordinate in self to log-pressure height (after duplicating self)
Return value: a GPhys
84 85 86 87 88 89 90 91 92 93 |
# File 'lib/numru/ganalysis/log_p.rb', line 84 def logp_coord_p2z(pdim=nil) pdim = GAnalysis::Met.find_prs_d(self) if !pdim p = self.coord(pdim) z = GAnalysis::LogP.p2z(p) ax = self.axis(pdim).copy ax.set_pos(z) ax.name = z.name grid = self.grid.copy.set_axis(pdim, ax) GPhys.new(grid,self.data) end |
#mouse_cut(dimx, dimy, num = 2, line_type = 1, line_index = 1) ⇒ Object
Makes a subset interactively by specifying a (poly-)line on the DCL viewport
ARGUMENTS
-
dimx {String] : name of number (0,1,..) of the dimension corresponding to the X coordinate in the current window of DCL
-
dimy {String] : name of number (0,1,..) of the dimension corresponding to the Y coordinate in the current window of DCL
-
num {Integer] : the number of points along the (poly-)line (2 or greater – if 2, a single line segment; if 3 or more, a poly-line)
RETURN VALUE
-
a GPhys
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 |
# File 'lib/numru/gphys/interpolate.rb', line 45 def mouse_cut(dimx, dimy, num=2, line_type=1, line_index=1) # < preparation > dimx = dim_index(dimx) dimy = dim_index(dimy) rundef = DCL.glpget("rundef") line = nil while(true) puts "\n*** Waiting for mouse click. Click #{num} points in the current viewport." line = DCLMouseLine.new(num) if line.ux.include?(rundef) puts "** The points specified include one(s) outside the U window. Do it again." else break end end line.draw(line_type, line_index) vx = line.vx vy = line.vy ux = line.ux uy = line.uy gpnew = dcl_fig_cut(dimx,dimy,ux,uy) [gpnew, line] end |
#mouse_cut_repeat ⇒ Object
Interpolation onto grid points specified by the previous call of GPhys#mouse_cut
169 170 171 172 173 174 175 176 |
# File 'lib/numru/gphys/interpolate.rb', line 169 def mouse_cut_repeat if @@interpo_previous_cutter.nil? raise("You must first use GPhys#mouse_cut. This method repeats it") end gpnew = interpolate(@@interpo_previous_cutter) gpnew = @@interpo_previous_modifier[gpnew] if @@interpo_previous_modifier gpnew end |
#phase_velocity(kdim, fdim, kconv, fconv, kf0_is_c0 = true, no_kfreorder = false) ⇒ Object
757 758 759 760 761 762 763 764 765 766 767 768 769 770 771 772 773 774 775 776 777 778 779 780 781 782 783 |
# File 'lib/numru/gphys/gphys_fft.rb', line 757 def phase_velocity(kdim,fdim,kconv,fconv,kf0_is_c0=true,no_kfreorder=false) kax = self.axis(kdim) fax = self.axis(fdim) kax.pos = kax.pos*kconv if kconv fax.pos = fax.pos*fconv if fconv cunits = fax.pos.units / kax.pos.units f = fax.pos.val k = kax.pos.val nk = k.length nf = f.length if no_kfreorder k[nk/2+1..-1] = -k[nk/2+1..-1][-1..0]+k[nk/2] f[nf/2+1..-1] = -f[nf/2+1..-1][-1..0]+f[nf/2] end f = -f cp = f.newdim(0) / k.newdim(1) #cp[kdim,fdim] jf0 = f.eq(0).where[0] # where f==0 jk0 = k.eq(0).where[0] # where k==0 if kf0_is_c0 cp[jk0,jf0] = 0.0 # treat k=f=0 as stationary (c=0) else cp[jk0,jf0] = 1.0/0.0 # not to count k=f=0 component at all (c=infty) end [cp, cunits] end |
#phase_velocity_binning(kdim, fdim, cbins, kconv = nil, fconv = nil) ⇒ Object
673 674 675 676 677 678 679 680 681 682 683 684 685 686 687 688 689 690 691 692 693 694 695 696 697 698 699 700 701 702 703 704 705 706 707 708 709 710 711 712 713 714 715 716 717 718 719 720 721 722 723 724 725 726 727 728 729 730 731 732 733 734 735 736 737 738 739 740 741 742 743 744 745 746 747 748 749 750 751 752 753 754 755 |
# File 'lib/numru/gphys/gphys_fft.rb', line 673 def phase_velocity_binning(kdim, fdim, cbins, kconv=nil, fconv=nil) # < process arguments > case cbins when Hash min = cbins["min"] ||raise(ArgumentError,"a Hash cbins must have 'min'") max = cbins["max"] ||raise(ArgumentError,"a Hash cbins must have 'max'") int = cbins["int"] ||raise(ArgumentError,"a Hash cbins must have 'int'") cbins = Array.new eps = int.abs*1e-6 # epsilon to deal with float steps (min.to_f..(max.to_f+eps)).step(int){|c| cbins.push(c)} cbins = NArray.to_na(cbins) when Array cbins = NArray.to_na(cbins) when NArray else raise ArgumentError, "cbins must be a Hash or Array or NArray" end kdim = dim_index(kdim) if kdim.is_a?(String) kdim += rank if kdim < 0 fdim = dim_index(fdim) if fdim.is_a?(String) fdim += rank if fdim < 0 # < sort along wavenumber/freuqency axis > pw = self.spect_zero_centering(kdim).spect_one_sided(fdim) # < process axes > cp, cunits = pw.phase_velocity(kdim,fdim,kconv,fconv,false) vcbins = VArray.new(cbins, {"units"=>cunits.to_s, "long_name"=>"phase velocity bounds"}, "cbounds") vccent = VArray.new( (cbins[0..-2] + cbins[1..-1])/2, {"units"=>cunits.to_s, "long_name"=>"phase velocity"}, "c") axc = Axis.new(true).set_cell(vccent, vcbins).set_pos_to_center axes = [axc] # the first dimension will be "c" gr = pw.grid (0...pw.rank).each do |d| if d!=kdim && d!=fdim axes.push(gr.axis(d)) end end newgrid = Grid.new(*axes) nk = pw.shape[kdim] nf = pw.shape[fdim] cp.reshape!(nk*nf) # < reorder input data > dimorder = (0...pw.rank).collect{|i| i} dimorder.delete(fdim) dimorder.unshift(fdim) dimorder.delete(kdim) dimorder.unshift(kdim) # --> [kdim, fdim, the other dims...] sh = pw.shape reshape = [nk*nf] (0...rank).each{|i| reshape.push(sh[i]) if i!=fdim && i!=kdim} pwv = pw.val.transpose(*dimorder).reshape(*reshape) # --> [ combined k&fdim, the other dims...] # < binning > shc = newgrid.shape pwc = NArray.new(pwv.typecode, *shc) # will have no missing data nc = axc.length for jc in 0...nc w = (cp.gt(cbins[jc]) & cp.lt(cbins[jc+1])).where pwc[jc,false] += pwv[w,false].sum(0) if w.length>0 w = (cp.eq(cbins[jc])).where pwc[jc,false] += pwv[w,false].sum(0)/2 if w.length>0 # half from bdry w = (cp.eq(cbins[jc+1])).where pwc[jc,false] += pwv[w,false].sum(0)/2 if w.length>0 # half from bdry end vpwc = VArray.new(pwc,pw.data,pw.name) gpwc = GPhys.new(newgrid,vpwc) gpwc end |
#phase_velocity_binning_iso_norml(kdim, fdim, cmin, cmax, cint, kconv = nil, fconv = nil) ⇒ Object
665 666 667 668 669 670 671 |
# File 'lib/numru/gphys/gphys_fft.rb', line 665 def phase_velocity_binning_iso_norml(kdim, fdim, cmin, cmax, cint, kconv=nil, fconv=nil) cbins = {"min"=>cmin,"max"=>cmax,"int"=>cint} pwc = phase_velocity_binning(kdim, fdim, cbins, kconv, fconv) fact = UNumeric[int, pwc.coord(0).units] pwc/fact end |
#phase_velocity_filter(xdim, tdim, cmin = nil, cmax = nil, xconv = nil, tconv = nil, remove_xtmean = false) ⇒ Object
628 629 630 631 632 633 634 635 636 637 638 639 640 641 642 643 644 645 646 647 648 649 650 651 652 653 654 655 656 657 658 659 660 661 662 663 |
# File 'lib/numru/gphys/gphys_fft.rb', line 628 def phase_velocity_filter(xdim, tdim, cmin=nil, cmax=nil, xconv=nil, tconv=nil, remove_xtmean=false) raise(ArgumentError,"need at least cmin or cmax") if !(cmin || cmax) xdim = dim_index(xdim) if xdim.is_a?(String) xdim += rank if xdim < 0 tdim = dim_index(tdim) if tdim.is_a?(String) tdim += rank if tdim < 0 fc = self.fft(nil,xdim,tdim) kdim = xdim fdim = tdim kconv = ( xconv ? 1.0/xconv : nil ) fconv = ( tconv ? 1.0/tconv : nil ) cp, = fc.phase_velocity(kdim,fdim,kconv,fconv,!remove_xtmean,true) fcv = fc.val nk = fc.shape[kdim] nf = fc.shape[fdim] sel = [true]*fc.rank for jf in 0...nf for jk in 0...nk c = cp[jk,jf] if ( cmin && c<cmin or cmax && c>cmax) sel[kdim]=jk sel[fdim]=jf fcv[*sel] = 0.0 end end end fc.replace_val(fcv) gp = fc.fft(true,xdim,tdim) gp = gp.real if (self.typecode <= NArray::FLOAT) GPhys.new(self.grid_copy, gp.data) #^ use the original grid, since units may have changed end |
#rawspect2powerspect(*dims) ⇒ Object
612 613 614 615 616 617 618 619 620 621 622 623 624 625 626 |
# File 'lib/numru/gphys/gphys_fft.rb', line 612 def rawspect2powerspect(*dims) # developpers memo: Needs Units conversion. factor = nil dims.each{|dim| ax = self.coord(dim) dwn = UNumeric.new( ((ax[-1].val - ax[0].val)/(ax.length - 1)).abs, ax.units ) if !factor factor = dwn**(-1) else factor = factor / dwn.to_f end } self * factor end |
#regrid(to) ⇒ Object
Interpolate to conform the grid to a target GPhys object
ARGUMENTS
-
to [GPhys] : the target gphys
RETURN VALUE
-
a GPhys
208 209 210 211 |
# File 'lib/numru/gphys/interpolate.rb', line 208 def regrid(to) coords = to.axnames.collect{|nm| to.coord(nm)} interpolate(*coords) end |
#running_mean(dim, len_or_wgt = nil, bc = BC_SIMPLE, nminvalid = 1) ⇒ Object
Running mean along a dimension (with optional weight specification).
ARGUMENTS
-
dim (Integer or String) : the dimension
-
len_or_wgt : If Integer, specifies the length; if 1D NArray, specifies the weight (e.g., NArray[1.0, 2.0, 1.0] for the 1-2-1 smooting)
-
bc (Integer; optional) : Speficy one of the folloing:
-
GPhys::BC_SIMPLE (default) : Averaging is trucated at the boundaries (the number of grid points used is reduced near the boundaries). The shape of the object is conserved.
-
GPhys::BC_CYCLIC : Cyclic boundary condition. Shape conserved.
-
GPhys::BC_TRIM : Grids near the boundaries are trimmed to secure the number of grid points to average. Shape not conserved; length along the dim is reduced by (len-1).
-
-
nminvalid (Integer; optional; defualt=1): This parameter is used only when the data have missing. Minimum number of grid points needed for averaging. Must be from 1 to len.
RETURN VALUE
-
a GPhys
REMARK AND LIMITATION
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 |
# File 'lib/numru/gphys/gphys_dim_op.rb', line 46 def running_mean(dim, len_or_wgt=nil, bc=BC_SIMPLE, nminvalid=1) #< process arguments > dim = dim_index(dim) # to handle String or negative specification case len_or_wgt when nil raise ArgumentError, "You need to specify the length (Integer) or the weight (1D NArray) as the 2nd argument" when Integer # len_or_wgt is a length len = len_or_wgt wgt = NArray.float(len).fill!(1.0) else # len_or_wgt is a weight wgt = len_or_wgt if (!wgt.respond_to?(:rank) || wgt.rank != 1) raise ArgumentError, "wgt: expect a 1D NArray(-like obj)" end len = wgt.length end #< calc running mean > vi = self.val if (vi.typecode > NArray::DFLOAT) raise("This method supports only real or integer data") end if vi.is_a?(NArrayMiss) vi, missval = nam2na_missval(vi) vo = c_running_mean(vi,dim,wgt,bc,missval,nminvalid) vo = NArrayMiss.to_nam(vo, vo.ne(missval) ) else vo = c_running_mean(vi,dim,wgt,bc) end #< grid > if (bc == BC_TRIM) fst = (len-1)/2 # if odd len/2, if even len/2-1 lst = -(len/2) - 1 grid = self.grid[ *([true]*dim + [fst..lst, false]) ] else grid = self.grid end #< result > vvo = VArray.new( vo, self.data, self.name ) # Inherit name & attrs GPhys.new( grid, vvo ) end |
#spect_one_sided(dim) ⇒ Object
601 602 603 604 605 606 607 608 609 610 |
# File 'lib/numru/gphys/gphys_fft.rb', line 601 def spect_one_sided(dim) dim = dim + self.rank if dim<0 len = self.shape[dim] b = self[ *([true]*dim + [0..len/2,false]) ] * 2 b[*([true]*dim + [0,false])] = b[*([true]*dim + [0,false])] / 2 if (self.shape[dim] % 2) == 0 # --> even number b[*([true]*dim + [-1,false])] = b[*([true]*dim + [-1,false])] / 2 end b end |
#spect_zero_centering(dim) ⇒ Object
587 588 589 590 591 592 593 594 595 596 597 598 599 |
# File 'lib/numru/gphys/gphys_fft.rb', line 587 def spect_zero_centering(dim) dim = dim + self.rank if dim<0 len = self.shape[dim] b = self[ *( [true]*dim + [[(len+1)/2..len-1,0..len/2],false] ) ].copy s1 = [true]*dim + [0, false] s2 = [true]*dim + [-1, false] if (len % 2) == 0 #--> even number b[*s1] = b[*s1]/2 # the ends are duplicated --> halved b[*s2] = b[*s1] end b.coord(dim)[0..len/2-1] = -b.coord(dim)[len/2+1..-1].val[-1..0] b end |
#threepoint_O2nd_deriv(*args) ⇒ Object
94 95 96 |
# File 'lib/numru/gphys/derivative.rb', line 94 def threepoint_O2nd_deriv(*args) Derivative::threepoint_O2nd_deriv(self,*args) end |