[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

[dennou-ruby:001966] Re: gpprint, gpview

竹広様
Cc: gtool4, dennou-ruby

塚原です.

> gphys をバージョンアップして無事 gpprint, gpview 動くようになりました. 
> dcl-c-gtk2 と dcl-c が conflict していて upgrade されてなかったようでした.

動いて良かったです. 

> で, 次々とすいませんが gpview について要望, 質問です. 
> 
>   ・アニメーション表示(gtview の -animate オプション)することは
>     まだできてないのでしょうか? 
> 
>   ・円筒座標, 地図投影しようと --itr オプションを使うと
> 
>      /usr/lib/ruby/1.6/numru/ggraph.rb:1270:in `grstrf': [GRSTRF] SIMFAC IS NOT DEFINED. (RuntimeError)
> 
>     となってしまいます. 地図投影パラメターを設定しないといけない
>     ようですが, どのように行うのでしょうか. 
> 

大変申しありません! 上記まだ対処できてないです. 特に円筒座標, 地図投影は
時間がかかりそうです. 来月頭までお待ちいただけないでしょうか?
# 大風呂敷広げておいて, 出来ないとは詐欺ですよね... すみません.

>     また, 2 次元図の座標軸の順序を入れ換えること(gtview での -exch オプション)
>     は難しいでしょうか. 

こちらは簡単でしたので実装してみました. --exch をつけて実行すると
gtview 同様に入換えます. 添付しましたのでよろしければ使ってみてください.
ついでに欠損値処理もするようにしましたので, 観測データも使えます.

# GGraph( GPhys 用描画モジュール ) が軸の入換えをサポートしているので
# とっても簡単でした. 素晴らしいです.

--------------------------------------
北海道大学院理学研究科 地球惑星科学専攻
地球流体力学研究室  M2  塚原大輔

email::daktu32@xxxxxxxxxxxxxxxxxxxx
--------------------------------------

#!/usr/bin/env ruby
##################################################
#=gpview
#
# quick viewer for the values of a variable specified by
# a gtool4-type URL.
#
#   2004/12/14  D Tsukahara && T Horinouti(parse_gturl)
#   2005/01/08  D Tsukahara (add option --exch and able to invalid value)
#==USAGE
#
#     % gpview url
#
# where the format of the url is
#
#     path@varname[,dimname=pos1[:pos2[:thinning_intv]][,dimname=...]]
#
# if you want to get more description, please see
#
#     http://www.gfd-dennou.org/arch/gtool4/gtool4-1.0/doc/gtview.htm 
#
#==EXAMPLES
#
#    % gpview data.nc@temp
#    % gpview data.nc@temp,lon=130:150,lat=0:90:2
#    % gpview --nocont data.nc@temp,lon=130:150,lat=0    
#    % gpview --noshade data.nc@temp,lon=130:150,lat=0    
#    % gpview --mean lon data.nc@temp,lon=130:150,lat=0
#    % gpview --exch data.nc@temp,lon=130:150,lat=0
#    
##################################################

require "getopts"        # for option_parse
require "numru/ggraph"   # ggraph library
include NumRu

#####################################################
$VIEWPORT = [0.15,0.85,0.23,0.58]
URLfmt = "path@varname[,dimname=pos1[:pos2[:thinning_intv]][,dimname=...]]"

## -- parse gturl --

def parse_gturl(gturl)
  if /(.*)@(.*)/ =~ gturl
    file = $1
    var = $2
  else
    raise "invalid URL: '@' between path & variable is not found" +
           "URL format: " + URLfmt
  end
  if /,/ =~ var
    slice = Hash.new
    thinning = Hash.new
    var_descr = var.split(/,/)
    var = var_descr.shift
    var_descr.each do |s|
      if /(.*)=(.*)/ =~ s
        dimname = $1
        subset = $2
        case subset
        when /(.*):(.*):(.*)/
          slice[dimname] = ($1.to_f)..($2.to_f)
          thinning[dimname] = {0..-1,$3.to_i}
        when /(.*):(.*)/
          slice[dimname] = ($1.to_f)..($2.to_f)
        else
          slice[dimname] = subset.to_f
        end
      else
        raise "invalid URL: variable subset specification error\n\n"
           "URL format: " + URLfmt
      end
    end
    thinning = nil if thinning.length == 0
  else
    slice = nil
    thinning = nil
  end
  [file, var, slice, thinning]
end

def parse_optanim(var)
  if /(.*),(.*)/ =~ var
    slice = Hash.new
    thinning = Hash.new
    var_descr = var.split(/,/)
    var = var_descr.shift
    var_descr.each do |s|
      if /(.*)=(.*)/ =~ s
        dimname = $1
        subset = $2
        case subset
        when /(.*):(.*):(.*)/
          slice[dimname] = ($1.to_f)..($2.to_f)
          thinning[dimname] = {0..-1,$3.to_i}
        when /(.*):(.*)/
          slice[dimname] = ($1.to_f)..($2.to_f)
        else
          slice[dimname] = subset.to_f
        end
      else
        raise "invalid range: variable subset specification error\n\n"
      end
    end
  end
  [var, slice, thinning]
end

def draw_map(gp, vp, style)
  l = gp.first2D.coord(0).val.min
  r = gp.first2D.coord(0).val.max
  b = gp.first2D.coord(1).val.min
  t = gp.first2D.coord(1).val.max
  DCL::grfig
  DCL::grstrn(10)
  DCL::grswnd(l,r,b,t)
  DCL::grsvpt(*vp)
  DCL::umpfit                      
  DCL::grstrf
  DCL::umpmap(style) 
end

def GGraph::annotate(str_ary)
  lnum = 0
  str_ary.each{ |str|lnum += 1 }
    charsize = 0.7 * DCL.uzpget('rsizec1')
  dvx = 0.01
  dvy = charsize*1.5
  raise TypeError,"Array expected" if ! str_ary.is_a?(Array)
  vxmin,vxmax,vymin,vymax = DCL.sgqvpt
  vx = 0.65
  vy = 0.045 + (lnum-1)*dvy
  str_ary.each{|str|
      DCL::sgtxzv(vx,vy,str,charsize,0,-1,1)
    vy -= dvy
  }
  nil
end

def draw(gp)
  # set missing value
  DCLExt.gl_set_params('lmiss'=>true,
                       'rmiss'=>gp.data.get_att('missing_value')[0]) if gp.data.get_att('missing_value')

  # fontsize
  DCL.sgpset('lcntl', false) ; DCL.uzfact(0.7)
  DCL.sgpset('lfull', true)  ; DCL.sgpset('lfprop',true)
  DCL.uscset('cyspos', 'B' )              # move unit y axis
  
  # viewport size
  GGraph.set_fig('viewport'=>$VIEWPORT)
  GGraph.set_fig('itr'=>($OPT_itr))
  GGraph.set_fig("xrev"=>"units:mb,units:hPa,units:millibar",  
                  "yrev"=>"units:mb,units:hPa,units:millibar") 
  
  # set options
  if $OPT_interval 
    GGraph.set_linear_contour_options('int'=>$OPT_interval)
    GGraph.set_linear_tone_options('int'=>$OPT_interval)
  end
  
  
  if ($OPT_line || gp.data.val.dim == 1) 
    draw_flag = "line"
  elsif (!$OPT_line && gp.data.val.dim >= 2) && !$OPT_noshade && $OPT_nocont 
    draw_flag = "nocont"
  elsif (!$OPT_line && gp.data.val.dim >= 2) && $OPT_noshade && !$OPT_nocont 
    draw_flag = "noshade"
  elsif (!$OPT_line && gp.data.val.dim >= 2) && !$OPT_noshade && !$OPT_nocont
    draw_flag = "full"
  end
  

  case draw_flag
  when "line"
    GGraph.line(gp, 
                true,
                "title"=>$OPT_title,
                "type"=>$OPT_type,
                "exchange"=>$OPT_exch
                )
  when "full"
    GGraph.tone(gp,
                true, 
                "title"=>$OPT_title,
                "transpose"=>$OPT_exch
                )
    GGraph.contour(gp, false, "transpose"=>$OPT_exch)
  when "nocont"
    GGraph.tone(gp,
                true, 
                "title"=>$OPT_title,
                "transpose"=>$OPT_exch
                )
  when "noshade"
    mj = DCL.udpget('indxmj')
    mn = DCL.udpget('indxmn')
    GGraph.contour(gp, 
                   true, 
                   "title"=>$OPT_title,
                   "label" =>true,
		   "transpose"=>$OPT_exch
                   )
  end
  
  if  draw_flag == ("full" || "nocont")
    #    color_bar#(if gphys colorbar options )
  end
  
  draw_map(gp, $VIEWPORT, $OPT_map) if $OPT_map
  
  
  DCL.ueitlv  
end

def draw_info(gp, file, var)
  DCL::sgtxzv($VIEWPORT[1]-0.07,$VIEWPORT[3]+0.02,"(units:#{gp.data.units.to_s})",
              0.7*DCL.uzpget('rsizec2'),0,0,3)    
  DCL::sgtxzv($VIEWPORT[1]-0.61,
              $VIEWPORT[3]+0.02,
              "[#{File::basename(file)}@#{var}]",
              0.7*DCL.uzpget('rsizec2'),0,0,3)  
end

def each_along_dims(gphys, loopdims, start)
  if !loopdims.is_a?(Array)
    loopdims = [loopdims]  # put in an Array (if a single Integer/String)
  end
  if loopdims.length == 0
    raise ArgumentError, "No loop dimension is specified "
  end
  
  
  loopdimids = Array.new
  loopdimnames = Array.new
  loopdims.each{|d|
    case d
    when Integer
      if d < 0
        d += gp.rank
	    end
      loopdimids.push( d )
      loopdimnames.push( gp.axis(d).name )
    when String
      loopdimids.push( gp.dim_index(d) )
      loopdimnames.push( d )
    else
      raise ArgumentError,"loopdims must consist of Integer and/or String"
    end
  }
  
  sh = Array.new
  len = 1
  loopdimids.each{|i|
    sh.push(gp.shape[i])
    len *= gp.shape[i]
  }
  
  cs = [1]
  (1...sh.length).each{|i| cs[i] = sh[i-1]*cs[i-1]}
  idx_hash = Hash.new
  for i in 0...len do
    loopdimnames.each_with_index{|d,j| 
      idx_hash[d] = ((i/cs[j])%sh[j])..((i/cs[j])%sh[j]) # rank preserved
    }
    subs = gphyses.collect{|g| g[idx_hash] }
    results = yield(*subs)
	  if !results.is_a?(Array)
	    raise "The return value of the block must be an Array of GPhys" 
	  end
    if i == 0
      results_whole = Array.new
      for j in 0...results.length
        rs = results[j]
        grid = rs.grid_copy
        loopdimnames.each{|nm|
          # replaces with original axes (full length)
          if !grid.axnames.include?( nm )
            raise "Dimension '#{nm}' has been eliminated. "+
                                                            "You must keep all loop dimensions." 
          end
          grid.set_axis(nm,gphyses[j].axis(nm))
        }
      end
    end
    idx = grid.axnames.collect{|k| idx_hash[k] || true}
    for j in 0...results.length
      rs = results[j]
      results_whole[j][*idx] = rs.data
    end
  end
  return results_whole
end


#####################################################
###++++++           Main Routine            ++++++###

unless getopts(
               "", 
               "wsn:", 
               "mean:", 
               "itr:1", 
               "title:", 
               "nocont", 
               "noshade", 
               "line", 
               "animate:", 
               "map:", 
               "interval:", 
	       "exch"
               )
  print "#{$0}:illegal options.\n"
  exit 1
end


# open netcdf var
gturl = ARGV[0]
file, var, slice,thinning = parse_gturl(gturl)
gp = GPhys::IO.open(file,var)
gp = gp.cut(slice) if slice
gp = gp[thinning]  if thinning

# mean along any axis
if ($OPT_mean)
  mean = ($OPT_mean).to_s.split(/\s*,\s*/)
  mean.each{|mn|
    gp = gp.mean(gp.axnames.index(mn))
  }
end

# draw figure
DCL.gropn($OPT_wsn||1)
#if $OPT_animate
#  each_along_dims(gp, 'lon', 0){|gp|
#    draw(gp)
#    draw_info(gp, file, var)
#  }
#else
  draw(gp)
  draw_info(gp, file, var)
#end
DCL.grcls