About this Blog

This Blog has English posts and Japanese posts. About Mac, iOS, Objective-C, and so on.

2013年4月25日木曜日

[ruby-opengl] 回転体 Advanced(回転角度によって半径を変えられる)

Z軸周りの回転体ならば、普通は、r=f(z)という関数を考えてそれをぐるっと回転させますが、
このメソッドでは、回転角度に応じて関数f(z)を変えて回転体を描くことができます。

例えば、theta=0のときはf(z)が定数で、theta=PIのときはf(z)は2次関数であるようにすれば、人のふくらはぎのような形も簡単に。



コードはこちら。
gl_cocoon.rb
require "opengl"

module GLUT
  def GLUT.Cocoon(height, slice, stack)
    return unless block_given?
    zStep = height/slice.to_f
    tStep = Math::PI/stack
    0.step(height-zStep+0.0001, zStep){|z|
       GL.Begin(GL::QUAD_STRIP)
      (-Math::PI).step(Math::PI+0.0001, tStep){|theta|
         r1 = yield(z, theta, tStep)
         r2 = yield(z+zStep, theta+tStep, tStep)
         phi = atan(zStep/(r2-r1).abs)
         GL.Normal(cos(theta)*sin(phi), sin(theta)*sin(phi), cos(phi))
         GL.Vertex(r1*cos(theta), r1*sin(theta), z)
         GL.Vertex(r2*cos(theta+tStep), r2*sin(theta+tStep), z+zStep)
      }
      GL.End
    }
  end


使い方
Z軸周りに回転します。
heightが高さ(Z軸方向の長さ)
sliceが高さ方向の刻み回数
stackが周方向の刻み回数
です。
渡すブロック内で、回転半径を返すようにしてください。

サンプル1
GLUT.Cocoon(6.0, 20, 20){|z,theta,tStep|
  # r = -a*z*(z-8) です。
  # 回転角が PI/3から PI*2/3の時はa=0.07、
  # それ以外の時は、a=0.1です。
  # 新宿のコクーンタワーのような形になります。
  case theta
  when Math::PI/3 .. 2*Math::PI/3
    -0.07*z*(z-8)
  else
    -0.1*z*(z-8)
  end
}



サンプル2 人の脛のような形も簡単に。
GLUT.Cocoon(6.0, 15, 15){|z,theta,tStep|
  # thetaは-PIからPIまで取るので、tは1〜0〜1と連続的に変化します。
  t = (theta/Math::PI).abs
  # tに1.0(定数)をかけ、1-tに-0.1*z^2 + 0.5*z + 1をかけることで、
  # f(z)は直線から2次曲線へと連続的に変化することになります。
  1.0*t + (1-t)*(-0.1*z*z + 0.5*z+1)
}


視点や光の設定などは省いています。このままでは表示されないので注意。

2013年4月23日火曜日

1つのベクトルから基本ベクトルを作る計算方法

なんかもうタイトルがわけわからないですが、
あるベクトルaが与えられた時に、そのベクトルとZ軸(またはX軸、Y軸)が対応しているような直交座標系(を表す単位ベクトル3つ)を求めるには?というお話です。

例によってRubyのプログラムで計算していきます。

require "matrix"
class Vector
  # 準備
  # Vectorクラスに、外積を求めるメソッドを追加します。
  # http://d.hatena.ne.jp/tanku/20100318/1268930228
  # よりコピペさせてもらいました。
  def outer_product(o)
    raise 'size#{self.size}' unless size == 3
    Vector[
      self[1] * o[2] - self[2] * o[1],
      self[2] * o[0] - self[0] * o[2],
      self[0] * o[1] - self[1] * o[0],
    ]
  end
  

  # 基本ベクトルを生成するメソッド。
  # ここでは、selfがZ軸に一致するようにしています。
  def generate_unit
    raise 'vector size must be 3.' unless size == 3
    # まず、selfと方向が一致する単位ベクトル w を作ります。
    w = self/self.r

    # 次に、wと直交するベクトル u を求めます。
    # 内積=0を満たしていればいいので、無数にあります。
    # ゼロ除算が起こらないような方法を選びました。
    # w[0]*w[1] - w[0]*w[1] + w[1]*w[2] - w[1]*w[2] + w[2]*w[0] - w[2]*w[0] = 0
    # を変形して、
    # w[0]*(w[1]-w[2]) + w[1]*(w[2]-w[0]) + w[2]*(w[0]-w[1]) = 0
    # としています。
    u0 = Vector[w[1]-w[2],  w[2]-w[0], w[0]-w[1]]
    u = u0/u0.r

    # wを軸として、uを90度回転したものが、残り1つのベクトル v です。
    # v = (w・u)*w + (u-(w・u)*w)*cos(90度) - (u✕w)*sin(90度)
    # 計算の理屈については、
    # http://hooktail.sub.jp/mathInPhys/vectorRot/
    # を見てください。
    v = w.inner_product(u)*w + u.outer_product(w)

    # 作った3つのベクトルを返します。
    return [u, v, w]
  end

参考にしたページのリンクを貼っておきます。
 3次元ベクトルメモ - 好き勝手に・げーあにん?
 ベクトルの回転 [物理のかぎしっぽ]
 class Vector

2013年4月21日日曜日

ruby-openglでサーモグラフィー風の色の生成

いろいろ計算方法はあるみたいだけど、このページの計算方法を試してみました。
 サーモグラフィー作成のためのRGB値の指定方法 - Ryota’s Research Diary


先に実行結果を


0.0から1.0までの値をとって、RGBを算出してカラーを設定します。

コードは以下。

gl_color.rb
module GL
  def GL.ThermoColor(t)
    case t
      when 0.0 .. 0.2
        r=0; g=0; b=5*t
      when 0.2 .. 0.4
        r=0; g=5*t-1; b=1
      when 0.4 .. 0.6
        r=0; g=1; b=3-5*t
      when 0.6 .. 0.8
        r=5*t-3; g=1; b=0
      when 0.8 .. 1.0
        r=1; g=5-5*t; b=0
      else 
        r=1; g=1; b=1
    end
    GL.Color(r,g,b)
  end
end
if $0 == __FILE__
  HEIGHT = 400
  WIDTH = 200
  def display
    GL.Clear(GL::COLOR_BUFFER_BIT)
    hStep = 2.0/100
    steps = 400
    steps.times{|t|
      h = 2*t.to_f/steps - 1.0
      GL.ThermoColor(t.to_f/steps)
      GL.Begin(GL::POLYGON)
      GL.Vertex(-1, h)
      GL.Vertex(1, h)
      GL.Vertex(1, h+hStep)
      GL.Vertex(-1, h+hStep)
      GL.End
    }
    GL.Flush()
  end
  GLUT.Init()
  GLUT.InitWindowSize(WIDTH, HEIGHT)
  GLUT.InitDisplayMode(GLUT::RGB)
  GLUT.CreateWindow("Thermo Color Sample")
  GLUT.DisplayFunc(method(:display).to_proc)
  GL.ClearColor(1.0, 1.0, 1.0, 1.0)  
  GLUT.MainLoop();
end

2013年4月16日火曜日

ruby-openglで楕円体を描いてみた

参考
http://www.gamedev.net/topic/126624-generating-an-ellipsoid-in-opengl/

2013/04/23: 修正しました。
・不要なインクルードを削除しました。
・GLU.NewQuadric()は必要なかったので、削除しました。
・GLUTの拡張ということにしました。(描画自体は、OpenGLのみで大丈夫です)

使い方
require "./gl_ellipsoid.rb"
# 中心は(0,0,0)
# 第2引数: X軸方向の半径
# 第3引数: Y軸方向の半径
# 第4引数: Z軸方向の半径
# 第5,第6引数: 描画の細かさ
GL.Ellipsoid(1, 2, 3, 20, 20)

gl_ellipsoid.rb
require "opengl"
require "glut"

module GLUT
  extend Math
  def GLUT.Ellipsoid(x,y,z,slice,stack)
    tStep = Math::PI/slice.to_f
    sStep = Math::PI/stack.to_f
#   ワイヤーフレームで描きたい時は、コメントを外してください。
#   GL.PolygonMode(GL::FRONT_AND_BACK, GL::LINE)
    (-Math::PI/2).step(Math::PI/2+0.0001, tStep){|t|
      GL.Begin(GL::TRIANGLE_STRIP)
      (-Math::PI).step(Math::PI+0.0001, sStep){|s|
        GL.Normal(cos(t)*cos(s), cos(t)*sin(s), sin(t))
        GL.Vertex(x*cos(t)*cos(s), y*cos(t)*sin(s), z*sin(t))
        GL.Vertex(x*cos(t+tStep)*cos(s), y*cos(t+tStep)*sin(s), z*sin(t+tStep))
      }
      GL.End
    }
  end
end


2013年4月14日日曜日

[ruby-opengl]Ubuntu12.10にインストール

3次元に挑戦!!ということで、ruby-openglを始めようと思います。

まったくRuby環境がインストールされていないなら、
$sudo apt-get install freeglut3 freeglut3-dbg freeglut3-dev
$sudo apt-get install ruby1.9.3
$sudo gem install opengl --pre
これで完了です。
--pre をつけるのは、このページのアドバイスに従いました。
Install ruby-opengl 0.60.1 on Ubuntu? ...in order to add 'gl', 'glu', and 'glut' extensions - Stack Overflow

以前は、gemのmkrfに不備があったようで、いろいろと面倒なことがあったらしいのですが、
(参考: Ubuntu 12.04.1 と ruby-1.9.3p283 でも ruby-opengl - Tociyuki::Diary
現在(2013/04/13)は、gem一発でいけました。


注意点
libopengl-ruby*がインストールされていると、干渉してうまく動かないようです。
Ubuntuソフトウエアセンターで、opengl rubyで検索して、インストールされていたら削除しておきましょう。


ちなみに、もうruby1.8が入っていて、rvmを使わずに1.9.3もインストールしたいなら、
$sudo apt-get install freeglut3 freeglut3-dbg freeglut3-dev
$sudo apt-get install ruby1.9.3
$sudo gem1.9.3 install opengl --pre

2013年4月12日金曜日

Rubyで画像の二値化

自習時のメモです。

参考
画像処理について
:: ロボティクス入門(共立出版株式会社)
BMP画像の形式について
:: http://www.kk.iij4u.or.jp/~kondo/bmp/
gnuplotについて
:: http://www.wakayama-u.ac.jp/~tokoi/lecture/shori1/2009/09.html

使い方
$ruby binarize.rb ファイル名
としてください。 入力は、グレースケールのWindowsビットマップ画像
二値化されたpbmファイルと、クラス間分散とクラス内分散の分散比の数値データを出力します。
数値データは、gnuplotでプロットしてください。

出力サンプル
元画像、二値化した画像、閾値を変化させたときの分散比、という順番です。

りんごの画像




机の上に置かれた、プラスチックのケース





コード
binarize.rb
if ARGV.size==0
  puts "ERROR!: input file name is missing."
  exit
end

$original_file_name = ARGV[0]

def output_pbm(data)
  height = data.size
  width = data[0].size
  File.open($original_file_name+".pbm", "w"){|io|
    io.puts "P1"
    io.puts "# CREATOR: Ruby PNM Filter Version 1.1"
    io.puts width.to_s + " " + height.to_s
    data.reverse.each{|e1|
      io.puts e1.join(" ")
    }
  }
  puts "Created #{$original_file_name+".pbm"}."
end
def get_bmp_data
  io = open($original_file_name, "rb")
  io.seek(2, IO::SEEK_CUR)
  print "file size (byte)="; puts filesize = io.read(4).unpack("i")[0]
  io.seek(4, IO::SEEK_CUR)
  print "offset="; puts offset = io.read(4).unpack("i")[0]
  io.seek(4, IO::SEEK_CUR)
  print "width="; p width = io.read(4).unpack("i")[0]
  print "height="; p height = io.read(4).unpack("i")[0]
  io.seek(2, IO::SEEK_CUR)
  print "Bits per pixel="; p bpp = io.read(2).unpack("s")[0]
  print "Compression type="; p ctype = io.read(4).unpack("i")[0]
  puts "(0:No compression 1:RunLength 8 bits/pixel 2:RunLength 4 bits/pixel 3:Bitfields)"
  io.seek(offset, IO::SEEK_SET)
  data = Array.new(height)
  height.times{|j|
    data[j] = Array.new(width)
    width.times{|i|
      value = io.read(1).unpack("C")[0]
      data[j][i] = value
    }
    io.seek(width%4==0 ? 0 : 4-width%4, IO::SEEK_CUR)
  }
  io.close
  return data
end
def binarize(data, threshold=128)
  return data.map{|e1| e1.map{|e2| e2 max_variance_ratio
      max_variance_ratio = v
      slashold  = s
    end
  }
}
puts "Created #{$original_file_name+".plt"} for gnuplot."
print "max_variance_ratio:"; puts max_variance_ratio
output_pbm(binarize(data, slashold))
print "time(sec): "; puts Time.now - start_time

2013年4月7日日曜日

[ソーシャルタギング]独自の共起尺度を考えてみた

きっかけ

ソーシャルブックマーク考2 タグの構造について より引用

上位概念Aと下位概念Bという二つの言葉がある場合、

①:「上位概念Aを含み、下位概念Bを含まない」文章( A && ! B )
②:「上位概念Aを含み、かつ下位概念Bを含む」文章( A & B)
③:「上位概念Aを含まず、下位概念Bを含む」文章( !A && B)

の3つの文章の数の比は、①>②>③が成り立つ。逆をいえば①>②>③が成り立っていれば、どちらが上位語かを判定できるわけだ。

上のサイトの仮説にもとづくと、
文章(webページ)の集合全体Uがあって、Aを含む文章、Bを含む文章の関係が図のようになっているときは、
Aが上位概念、Bが下位概念である、と言えます。
この状態を便宜的に「状態I」とします。


AとBの大小を固定しておいた場合、①>②>③という大小関係の他に、
①>③>②という場合も考えられます。
つまり、下の図のようになっているときです。この場合、AとBの関連は薄いと言えます。
この状態を便宜的に「状態II」とします。


発展

しかし、実際問題としては「①>②>③」と「①>③>②」をすっぱりと判定できるものではありません。
そこで、「A & B」を、A、Bそれぞれで割ることで正規化し、座標上にプロットしてみます。 (しばらく、NA>NBの場合を考えます。)

,


すると、状態Iは 


各辺をで割って、


つまり、 となり、

状態Iのときは上の網掛けの正方形部分にプロットされます。


同様に、状態IIのときは
が、
となり、
となります。

状態IIは上の網掛けの三角形部分にプロットされます。


そうなると気になるのが、

上の図の網掛け部分にプロットされるような場合です。

,
ベン図で表すなら、このような、AとBが重なってしまうような状況。

おそらく、文章を書く上での「揺れ」とみなされるような場合であると考えられます。
つまり、AとBは同意語と言えます。


(いったん)まとめ

2つのタグA,Bの関係を調べるのに、
,
を考えて座標上にプロットします。

a.にある時・・・AとBの関連は薄い
b.にある時・・・AはBの上位概念である
c.にある時・・・AとBは同意語


さらに考察

が、境界ですっぱりと分けてしまうのは多少乱暴な気もします。
「傾向が強い」「弱い」「可能性が高い」「低い」としておいたほうがよさそうです。
そこで、原点(0,0)からプロットした点(x,y)への距離Rと、
y=xと(0,0)-(x,y)のなす角度を考えてみます。


仮定

Rが大きいほど、AとBの関連度が高い。0〜R〜ルート2。
角度(θ)が+に大きいほど、AはBよりも上位の概念である可能性が高い。-45度〜θ〜45度

まとめ

タグの上位下位を決定できるような仮説を立てることができました。
実際にはてなブックマークのデータを使って実験したものがあるので、そのうち公開したいと思います。

[ソーシャルタギング]共起性尺度のまとめ

以前、taglaut.comというドメイン(現在はありません)で行なっていた、ソーシャルタギングのメモを発見したので公開しておきます。

共起の例と応用

ある人が、通販サイトからAという商品とBという商品を買った

多くの人が商品Aと商品Bを買っているなら、共起性尺度が高い、と言えます。
そこで、Aを買った人にBをおすすめ商品として紹介すれば、買ってもらえる可能性は高いです。
(amazon.comをよく使う人ならピンとくるでしょう)
出てくる結果はある意味で当たり前のことなので、それほどすごいとは感じないかもしれません。
例えば、釣竿を買う人なら釣り糸や餌・ルアー、釣り用ベストなどをおすすめすればいいというのは常識として分かります。
これを購買データの解析で自動で行えるというのがすごいところです。
また、誰も思いつかないような傾向が出てくるかもしれません。


ソーシャルブックマークで、あるページに対して、AというタグとBというタグが付けられた。

タグAとタグBの共起性尺度が高いなら、あるページにタグAがついた時点で、タグBもそのページにふさわしいだろうということが予想出来ます。


いろいろな共起性尺度

ここでは、ソーシャルブックマークでのタグ付けを念頭においています。
UがWebページの全体、
AがタグAがつけられたページの集合、NAがその個数、
BがタグBがつけられたページの集合、NBがその個数
です。
(同じページに同じタグが複数回付くこともありますが、ここでは考えていません。)

Jaccard係数

「AとBが共起した回数」割る「AとBの和集合」で求められます。

Simpson係数

「AとBが共起した回数」割る「A,Bの回数のうちの小さい方」で求められます。

Dice係数

「AとBが共起した回数」割る「A,Bの回数の平均」で求められます。


参考

[1]統計ソフトRのブログ 共起性尺度
[2]WEB+DB PRESS Vol.57 アルゴリズム実践教室【第2回】レコメンドエンジン開発に挑戦―関連記事を導き出すしくみを知る 伊藤直也
[3]類似度と距離 - CatTail Wiki*
[+]数式の表現にはCODECOGSのサービスを使っています。