vtkSphereSourceは便利だが,勉強にならないので,球をセルから構成してみる.
定数
package require vtk set nlon 100 set nlat [expr {$nlon/2 + 1}] set pi [expr {acos(-1.)}] set pi2 [expr {$pi * 2.}] set npts [expr {$nlon * $nlat}] #puts stdout "pi=$pi pi2=$pi2 npts=$npts"
経度φと余緯度θからデカルト座標に変換するコマンドを作る. VTKでは,経度がθ余緯度がφのようだが気にしない.
proc phitheta2xyz args { set phi [lindex $args 0] set theta [lindex $args 1] set xyz {} set sin_theta [expr {sin($theta)}] lappend xyz [expr {cos($phi) * $sin_theta}] lappend xyz [expr {sin($phi) * $sin_theta}] lappend xyz [expr {cos($theta)}] }
経度,余緯度を生成.
set lon {} set dlon [expr {$pi/$nlon}] for {set i 0} {$i<$nlon} {incr i} { lappend lon [expr {$pi2 * $i / $nlon}] } #puts stdout $lon set lat {} set dlat [expr {$pi/($nlat-1)}] for {set j 0} {$j<$nlat} {incr j} { lappend lat [expr {$dlat * $j}] } #puts stdout $lat
座標点のデカルト座標を生成してvtkPointsに登録.
vtkPoints pts pts SetNumberOfPoints $npts set i 0 foreach theta $lat { foreach phi $lon { # puts stdout "$phi $theta" # puts [phitheta2xyz $phi $theta] eval pts SetPoint [concat $i [phitheta2xyz $phi $theta]] incr i } }
4点からセルを作る.経度方向に周期境界.$nlon-1番目の隣は,0番目に戻る.
vtkCellArray clls set nlonm1 [expr {$nlon - 1}] for {set j 0} {$j<$nlat-1} {incr j} { for {set i 0} {$i<$nlon} {incr i} { set k0 [expr {$nlon*$j + $i}] if {$i!=$nlonm1} { set k1 [expr {$k0 + 1}] } else { set k1 [expr {$k0 - $nlonm1}] } set k2 [expr {$k1 + $nlon}] set k3 [expr {$k0 + $nlon}] clls InsertNextCell 4 clls InsertCellPoint $k0 clls InsertCellPoint $k1 clls InsertCellPoint $k2 clls InsertCellPoint $k3 } }
点とセルをvtkPolyDataに登録.
vtkPolyData sfc sfc SetPoints pts sfc SetPolys clls pts Delete clls Delete
ファイルに保存して完了.
vtkPolyDataWriter writer writer SetInput sfc writer SetFileName sphere2.vtk writer Write writer Delete sfc Delete