Posted to tcl by ahype18 at Thu Jan 08 22:08:32 GMT 2015view raw
- # Routine that extracts a VOI based on the sections and integrates a volume based on
- # a the extracted volume, a cutoff and the elevations planes
- #
- set showvol 1
- set savestrgrd 0
- toplevel .contour
- frame .contour.buttons -borderwidth 3 -relief ridge
- pack .contour.buttons -side top -padx 8 -pady 5 -ipadx 3 -ipady 1
- frame .contour.buttons.f1
- label .contour.buttons.l1 -text "Thresh Upper:" -width $lwidth -anchor w
- entry .contour.buttons.r1 -textvariable volcontvalU -width $rwidth
- pack .contour.buttons.f1
- pack .contour.buttons.l1 .contour.buttons.r1 -in .contour.buttons.f1 -side left -padx $xpad -pady $ypad
- frame .contour.buttons.f1a
- label .contour.buttons.l1a -text "Thresh Lower:" -width $lwidth -anchor w
- entry .contour.buttons.r1a -textvariable volcontvalL -width $rwidth
- pack .contour.buttons.f1a
- pack .contour.buttons.l1a .contour.buttons.r1a -in .contour.buttons.f1a -side left -padx $xpad -pady $ypad
- frame .contour.buttons.f1b
- label .contour.buttons.l1b -text "BVWirr:" -width $lwidth -anchor w
- entry .contour.buttons.r1b -textvariable bvwirr -width $rwidth
- pack .contour.buttons.f1b
- pack .contour.buttons.l1b .contour.buttons.r1b -in .contour.buttons.f1b -side left -padx $xpad -pady $ypad
- frame .contour.buttons.f1c
- label .contour.buttons.l1c -text "Bo:" -width $lwidth -anchor w
- entry .contour.buttons.r1c -textvariable shrinkage -width $rwidth
- pack .contour.buttons.f1c
- pack .contour.buttons.l1c .contour.buttons.r1c -in .contour.buttons.f1c -side left -padx $xpad -pady $ypad
- frame .contour.buttons.f2
- button .contour.buttons.f2.b1 -text "Contour" -width 10 -command contour
- pack .contour.buttons.f2
- pack .contour.buttons.f2.b1 -side left -padx 4 -pady 4
- frame .contour.save -borderwidth 3 -relief ridge
- pack .contour.save -side top -padx 8 -pady 5 -ipadx 3 -ipady 1
- frame .contour.save.f1
- label .contour.save.l1 -text "Save Poly:" -width $lwidth -anchor w
- entry .contour.save.r1 -textvariable volfile -width $rwidth
- pack .contour.save.f1
- pack .contour.save.l1 .contour.save.r1 -in .contour.save.f1 -side left -padx $xpad -pady $ypad
- frame .contour.save.f2
- button .contour.save.f2.b1 -text "Save" -width 10 -command volsave
- pack .contour.save.f2
- pack .contour.save.f2.b1 -side left -padx 4 -pady 4
- frame .contour.save.f1b
- label .contour.save.l1b -text "Save Map:" -width $lwidth -anchor w
- entry .contour.save.r1b -textvariable mapfile -width $rwidth
- pack .contour.save.f1b
- pack .contour.save.l1b .contour.save.r1b -in .contour.save.f1b -side left -padx $xpad -pady $ypad
- frame .contour.save.f2b
- button .contour.save.f2b.b1 -text "Save" -width 10 -command mapsave
- pack .contour.save.f2b
- pack .contour.save.f2b.b1 -side left -padx 4 -pady 4
- frame .contour.save.f1c
- label .contour.save.l1c -text "Save StrGrd:" -width $lwidth -anchor w
- entry .contour.save.r1c -textvariable strgrdfile -width $rwidth
- pack .contour.save.f1c
- pack .contour.save.l1c .contour.save.r1c -in .contour.save.f1c -side left -padx $xpad -pady $ypad
- frame .contour.save.f2c
- button .contour.save.f2c.b1 -text "Save" -width 10 -command strgrdsave
- pack .contour.save.f2c
- pack .contour.save.f2c.b1 -side left -padx 4 -pady 4
- radiobutton .contour.buttons.rb1 -variable threshtype -text "Higher" -value "1"
- pack .contour.buttons.rb1 -side left
- radiobutton .contour.buttons.rb2 -variable threshtype -text "Between" -value "0"
- pack .contour.buttons.rb2 -side left
- radiobutton .contour.buttons.rb3 -variable threshtype -text "Lower" -value "-1"
- pack .contour.buttons.rb3 -side left
- checkbutton .contour.buttons.cb1 -variable showvol -text "Show Vol" -command showvol
- pack .contour.buttons.cb1 -side top -anchor w
- text .contour.t -height 20 -width 60 -background white
- pack .contour.t
- set volcontvalU 0.0
- set volcontvalL 0.0
- set bvwirr 0.06
- set shrinkage 1.1
- set volrunflg 0
- set volfile "VolPoly.vtk"
- set mapfile "Map.dat"
- set strgrdfile "StrGrd-Extract.vtk"
- proc contour {} {
- global volcontvalU volcontvalL status dx dy gx1 gx2 gy1 gy2 z1 z2 elev1 elev2 volrunflg bvwirr shrinkage sumvol cellvol
- global grdindex sumnetthick grdsize volxdim volxmin volymin threshtype passcut sumpor ptscal sw sumhc zgriddim sumthick
- set status "Contouring" ; update
- .contour.t delete 1.0 end
- if {$volrunflg>0} {
- contGeom Delete
- p2c Delete
- thresh Delete
- connect Delete
- transConn Delete
- lut2 Delete
- geomFilt Delete
- connTriang Delete
- connSmooth Delete
- connNormals Delete
- moldMapper Delete
- ren1 RemoveActor moldActor
- moldActor Delete
- }
- vtkExtractGrid contGeom
- contGeom SetInput [reader GetOutput]
- contGeom SetVOI $gx1 $gx2 $gy1 $gy2 $z1 $z2
- contGeom Update
- set voldim [[contGeom GetOutput] GetDimensions]
- set volxdim [lindex $voldim 0]
- set volydim [lindex $voldim 1]
- set volzdim [lindex $voldim 2]
- .contour.t insert end "GX1= $gx1\n"
- .contour.t insert end "GX2= $gx2\n"
- .contour.t insert end "GY1= $gy1\n"
- .contour.t insert end "GY2= $gy2\n"
- .contour.t insert end "Zmin= $z1\n"
- .contour.t insert end "Zmax= $z2\n"
- .contour.t insert end "Elev1= $elev1\n"
- .contour.t insert end "Elev2= $elev2\n"
- .contour.t insert end "Xdim= $volxdim\n"
- .contour.t insert end "Ydim= $volydim\n"
- .contour.t insert end "Zdim= $volzdim\n"
- set ncells [expr $volxdim*$volydim*$volzdim]
- .contour.t insert end "nCells= $ncells\n"
- set passcut 0
- set sumpor 0.0
- set sumhc 0.0
- set avgpor 0.0
- set sumthick 0.0
- set sumvol 0.0
- set grdsize [expr $volxdim*$volydim]
- set repor 0.0
- for {set i4 0} {$i4 < $grdsize} {incr i4} {
- set sumnetthick($i4) 0.0
- }
- set volpt0 [[contGeom GetOutput] GetPoint 0]
- set volxmin [lindex $volpt0 0]
- set volymin [lindex $volpt0 1]
- for {set i5 0} {$i5 < $ncells} {incr i5} {
- set calck [expr $i5/$grdsize]
- set grdindex [expr $i5-$calck*$grdsize]
- set volpt [[contGeom GetOutput] GetPoint $i5]
- set volz [lindex $volpt 2]
- if {$calck < [expr $volzdim-1]} {
- set i5b [expr $i5 + $grdsize]
- } else {
- set i5b [expr $i5 - $grdsize]
- }
- set volptb [[contGeom GetOutput] GetPoint $i5b]
- set volzb [lindex $volptb 2]
- set zgriddim [expr abs($volz - $volzb)]
- set cellvol [expr $dx*$dy*$zgriddim]
- set ptscal [[[[contGeom GetOutput] GetPointData] GetScalars] GetValue $i5]
- if {$volz>=$elev1 && $volz<=$elev2} {
- if {$threshtype==1 && $ptscal>= $volcontvalU} {
- sumcells
- } elseif {$threshtype==0 && $ptscal>=$volcontvalL && $ptscal<=$volcontvalU} {
- sumcells
- } elseif {$threshtype==-1 && $ptscal<=$volcontvalL} {
- sumcells
- }
- } else {
- eval [[[[contGeom GetOutput] GetPointData] GetScalars] SetValue $i5 $repor]
- }
- }
- #set avgpor [expr $sumpor/$sumvol]
- #set porevolbbls [expr $sumpor*0.178108]
- #set hcvolbbls [expr $sumhc*0.178108]
- #.contour.t insert end "PassCont= $passcut\n"
- #.contour.t insert end "AvgPor= $avgpor\n"
- #.contour.t insert end "PoreVol= $sumpor\n"
- #.contour.t insert end "PoreVol(bbl)= $porevolbbls\n"
- #.contour.t insert end "HCVol= $sumhc\n"
- #.contour.t insert end "HCVol(bbl)= $hcvolbbls\n"
- update
- vtkPointDataToCellData p2c
- p2c SetInput [contGeom GetOutput]
- vtkThreshold thresh
- thresh SetInput [contGeom GetOutput]
- if {$threshtype==1} {thresh ThresholdByUpper $volcontvalU}
- if {$threshtype==0} {thresh ThresholdBetween $volcontvalL $volcontvalU}
- if {$threshtype==-1} {thresh ThresholdByLower $volcontvalL}
- vtkConnectivityFilter connect
- connect SetInput [thresh GetOutput]
- connect SetExtractionModeToAllRegions
- connect ColorRegionsOn
- connect Update
- set scalrange2 [[connect GetOutput] GetScalarRange]
- set scalmin2 [lindex $scalrange2 0]
- set scalmax2 [lindex $scalrange2 1]
- .debug.t insert end "Scalmin2= $scalmin2\n"
- .debug.t insert end "Scalmax2= $scalmax2\n" ; update
- vtkLookupTable lut2
- lut2 SetHueRange 0.667 0.0
- lut2 SetTableRange $scalmin2 $scalmax2
- lut2 Build
- vtkGeometryFilter geomFilt
- geomFilt SetInput [connect GetOutput]
- vtkTriangleFilter connTriang
- connTriang SetInput [geomFilt GetOutput]
- vtkSmoothPolyDataFilter connSmooth
- connSmooth SetInput [connTriang GetOutput]
- connSmooth BoundarySmoothingOn
- vtkTransformFilter transConn
- transConn SetInput [connSmooth GetOutput]
- transConn SetTransform trans1
- transConn Update
- vtkPolyDataNormals connNormals
- connNormals SetInput [transConn GetOutput]
- vtkPolyDataMapper moldMapper
- moldMapper SetInput [connNormals GetOutput]
- moldMapper ScalarVisibilityOn
- moldMapper SetLookupTable lut2
- moldMapper SetScalarRange $scalmin2 $scalmax2
- vtkActor moldActor
- moldActor SetMapper moldMapper
- set volrunflg 1
- showvol
- set status "Ready" ; update
- }
- proc sumcells {} {
- global passcut sumpor ptscal sw bvwirr sumhc shrinkage sumnetthick grdindex zgriddim sumthick sumvol cellvol
- incr passcut
- set sumthick [expr $sumthick+$zgriddim]
- set sumvol [expr $sumvol+$cellvol]
- set sumpor [expr $sumpor+$ptscal*$cellvol]
- if {$ptscal <= 0.0} {
- set sw 1.0
- } else {
- set sw [expr $bvwirr/$ptscal]
- }
- set sumhc [expr $sumhc+$cellvol*$ptscal*(1-$sw)/$shrinkage]
- set sumnetthick($grdindex) [expr $sumnetthick($grdindex)+$zgriddim] ; # feet meeting cutoff
- # set sumnetthick($grdindex) [expr $sumnetthick($grdindex)+$zgriddim*(1.0-$ptscal)] ; # net sand thIckness
- # set sumnetthick($grdindex) [expr $sumnetthick($grdindex)+$zgriddim*$ptscal] ; # porosity thickness
- }
- proc showvol {} {
- global showvol
- set status "Displaying Volume" ; update
- set existsVol [[ren1 GetActors] IsItemPresent moldActor]
- if {$showvol==0 && $existsVol>0} {ren1 RemoveActor moldActor}
- if {$showvol==1 && $existsVol<=0} {ren1 AddActor moldActor}
- renWin Render
- if {$showvol==0 && $existsVol>0} {
- ren1 RemoveActor moldActor
- } elseif {$showvol==1 && $existsVol<=0} {
- ren1 AddActor moldActor
- }
- }
- set volsaveflg 0
- proc volsave {} {
- global status volsaveflg volfile
- set status "Saving Vol Poly" ; update
- if {$volsaveflg>0} {writeVolPoly Delete}
- vtkPolyDataWriter writeVolPoly
- writeVolPoly SetInput [connSmooth GetOutput]
- writeVolPoly SetFileName data/$volfile
- writeVolPoly SetFileTypeToBinary
- writeVolPoly Update
- set volsaveflg 1
- set status "Ready" ; update
- }
- proc mapsave {} {
- global status sumnetthick mapfile grdsize volxdim dx dy volxmin volymin
- set status "Saving Map Grid" ; update
- set mapFileId [open $mapfile w]
- for {set i 0} {$i < $grdsize} {incr i} {
- # set calcgrdj [expr double($i/$volxdim)]
- # set grdj [expr double(int($calcgrdj)*$dy)+$volymin]
- # set calcgrdi [expr $calcgrdj*double($volxdim)]
- # set grdi [expr double(($i-$calcgrdi)*$dx)+$volxmin]
- set grdpt [[contGeom GetOutput] GetPoint $i]
- set grdx [lindex $grdpt 0]
- set grdy [lindex $grdpt 1]
- puts $mapFileId "$grdx $grdy $sumnetthick($i)"
- }
- close $mapFileId
- set status "Ready" ; update
- }
- set savestrgrdflg 0
- proc strgrdsave {} {
- global status strgrdfile savestrgrdflg gx1 gx2 gy1 gy2 z1 z2
- set status "Saving Structured Grid" ; update
- if {$savestrgrdflg==1} {
- writeStrgrd Delete
- }
- vtkExtractGrid contGeom
- contGeom SetInput [reader GetOutput]
- contGeom SetVOI $gx1 $gx2 $gy1 $gy2 $z1 $z2
- contGeom Update
- vtkStructuredGridWriter writeStrgrd
- writeStrgrd SetInput [contGeom GetOutput]
- writeStrgrd SetFileName $strgrdfile
- # writeStrgrd SetFileTypeToBinary
- writeStrgrd Update
- set savestrgrdflg 1
- contGeom Delete
- set status "Ready" ; update
- }