gis
gis
管理员
管理员
  • 注册日期2003-07-16
  • 发帖数15947
  • QQ554730525
  • 铜币25339枚
  • 威望15364点
  • 贡献值0点
  • 银元0个
  • GIS帝国居民
  • 帝国沙发管家
  • GIS帝国明星
  • GIS帝国铁杆
阅读:2400回复:7

不管了,贴个模块,大家看看先!

楼主#
更多 发布于:2003-07-28 14:23
'
' ESRI
' 3D Analyst Developer Sample -
' bas3DArea.bas

' Contains a method 'UpdatePolysWith3DSurfaceArea' for deriving the 3D surface area
' for polygons by using a TIN edit session

'
' May require references to:
' ESRI ArcScene Object Library
' ESRI ArcMap Object Library
' ESRI TIN Object Library
' ESRI Object Library
'
Option Explicit

'
' looping through each polygon in the designated layer:
' 1. clip and drape it on the TIN to create a 3D poly
' 2. for each part of the polygon (consider multipart geometries),
' 3. burn the geometry into the TIN using a unique tag value for the underlying triangles
' 4. extract the new underlying triangles under the original geometry
' 5. use these triangles as input for the GetPartialVolumeAndArea to calculate the area for the geometry
' 6. update the designated field of the polygon feature class with the 3D area
'
Public Sub UpdatePolysWith3DSurfaceArea(sTINLayer As String, sPolyLayer As String, sPolyFieldToUpdate As String)

  On Error GoTo UpdatePolysWith3DSurfaceArea_ERR
  
  Dim pTinAdv As ITinAdvanced
  Dim pTinLayer As ITinLayer
  Dim pFC As IFeatureCursor
  Dim pFeature As IFeature
  Dim pPoly As IPolygon
  Dim pTinEdit As ITinEdit
  Dim pTinFEdit As ITinFeatureEdit
  Dim pFilt As ITinFilter
  Dim pRow As IRow
  Dim pTinSurface As ITinSurface
  Dim nArea As Double, nProjArea As Double, nRef As Double
  Dim pEnumTriangles As IEnumTinTriangle
  Dim pWksEdit As IWorkspaceEdit
  Dim pWks As IWorkspace
  Dim pWksFactory As IWorkspaceFactory
  Dim nUpdateField As Long
  Dim nFeatures As Long
  
  ' access the TIN layer and the polygon layer:
  Set pTinLayer = bas3DArea.GetLayer(sTINLayer)
  If pTinLayer Is Nothing Then
    Debug.Assert 0
    MsgBox "The surface: " & sTINLayer & " was not found.", vbCritical, "UpdatePolysWith3DSurfaceArea"
    Exit Sub
  End If
    
  Set pFC = bas3DArea.GetFeatureCursorFromLayer(sPolyLayer)
  If pFC Is Nothing Then
    Debug.Assert 0
    MsgBox "A feature cursor could not be established from: " & sPolyLayer, vbCritical, "UpdatePolysWith3DSurfaceArea"
    Exit Sub
  End If
  
  ' reference the polygon's field to be updated:
  Dim pFLayer As IFeatureLayer
  Set pFLayer = bas3DArea.GetLayer(sPolyLayer)
  nUpdateField = pFLayer.FeatureClass.FindField(sPolyFieldToUpdate)
  If nUpdateField < 0 Then
    Debug.Assert 0
    MsgBox "The field '" & sPolyFieldToUpdate & "' was not found in: " & sPolyLayer & ".", vbCritical, "UpdatePolysWith3DSurfaceArea"
    Exit Sub
  End If
  
  ' reference necessary interfaces:
  Set pTinAdv = pTinLayer.Dataset
  Set pTinEdit = pTinAdv
  Set pTinFEdit = pTinAdv
  Set pTinSurface = pTinAdv
  
  ' use a reference plane of the z extent of the TIN since we will
  ' be getting area from selected triangles below it:
  nRef = pTinAdv.Extent.zmax

  ' looping through each polygon in the feature class, start temporary in-memory
  ' edit sessions which will:
  ' 1. clip and drape it on the TIN to create a 3D poly
  ' 2. for each part of the polygon (consider multipart geometries),
  ' 3. burn the geometry into the TIN using a unique tag value for the underlying triangles
  ' 4. extract the new underlying triangles under the original geometry
  ' 5. use these triangles as input for the GetPartialVolumeAndArea to calculate the area for the geometry
  ' 6. update the designated field of the polygon feature class with the 3D area
  Set pFeature = pFC.NextFeature
  If Not TypeOf pFeature.Shape Is IPolygon Then
    MsgBox "Error: " & sPolyLayer & " does not contain polygons.", vbCritical, "UpdatePolysWith3DSurfaceArea"
    Exit Sub
  End If
  
  Dim bUseFeature As Boolean
  
  ' loop for each polygon feature, calculating area values by converting the
  ' polys to 3D, 'burning' into the TIN during an edit session to establish new
  ' triangles, from which to query and use as input to the ITinSurface::GetPartialAreaAndVolume method:
  Do While Not pFeature Is Nothing
  
    bUseFeature = True
    If Not pTinFEdit.StartInMemoryEditing() Then
      MsgBox "An in-memory edit session could not be established for: " & sTINLayer & ".", vbCritical, "UpdatePolysWith3DSurfaceArea"
      Exit Sub
    End If
    
    Set pPoly = pFeature.Shape
    
    ' create a 3D polygon, clipped if necessary to the surface extent:
    Dim p3DPoly As IPolygon2
    pTinSurface.InterpolateShape pPoly, p3DPoly
    
    If p3DPoly Is Nothing Then
      bUseFeature = False
    End If
    
    ' Get polygon connected components - to handle multipart polygons; those with one or
    ' more external and/or internal rings. The TinFeatureEdit calls used later rely on
    ' simple connected components.
    Dim lExternalRings As Long
    lExternalRings = p3DPoly.ExteriorRingCount

    If (lExternalRings > 50) Then ' TODO - note the hardcode limit
      ' MsgBox "This command can not honor multipart geometries greater than 50 parts.", vbExclamation, "UpdatePolysWith3DSurfaceArea"
      bUseFeature = False
      Exit Sub
    End If
  
    If bUseFeature Then
      ' bUseFeature = true
 
      ' get the individual parts of the polygon-
      ' of course, there may be only 1 in the case of a simple polygon,
      ' but this way, inside a loop of each part, we will support non-simple polygons as well
      ' for this area calculation:
      Dim pConnectedComps(0 To 49) As IPolygon
      p3DPoly.GetConnectedComponents lExternalRings, pConnectedComps(0)
        
      Dim pInterpGC As IGeometryCollection
      Set pInterpGC = p3DPoly
    
      ' generate a unique tag value in which to tag all underlying triangles under the
      ' geometry we will temporarily add to the TIN, and create a 'seed' which
      ' will allow us to extract all triangles with the similiar tags-
      ' this is leading up to the creation of a triangle enumerator which
      ' will serve as the input to our ITinSurface:GetPartialAreaAndVolume method:
      Dim lTag As Long
      Dim pSeed As ITinFeatureSeed
      Dim pFilter As ITinValueFilter
      
      lTag = pTinAdv.GenerateUniqueTagValue(esriTinTriangle)
      Set pSeed = New TinTriangle
      pSeed.UseTagValue = True
      
      Set pFilter = New TinValueFilter
      pFilter.ActiveBound = esriTinUniqueValue
      pFilter.UniqueValue = lTag
    
      ' in a loop, add each polygon geometry to the TIN, tagging the triangles formed.
      ' then extract all those triangles with the common tag, and use as the input
      ' to the surface area method.
      ' this should be the equivalent of the 3D area of the 2D polygon used to
      ' tag these triangles:
      Dim i As Long
      nArea = 0
      For i = 0 To (lExternalRings - 1)
          
        Dim pSimplePoly As IPolygon
        Set pSimplePoly = pConnectedComps(i)
            
        pTinFEdit.AddPolygonZ pSimplePoly, esriTinHardEdge, lTag, 0, 0, pSeed
    
        Dim pTinPoly As ITinPolygon
        Set pTinPoly = pTinAdv.ExtractPolygon(pSeed, pFilter, False)
        
        Set pEnumTriangles = pTinPoly.AsTriangles
      
        Dim nAreaTemp As Double
        ' call the ITinSurface::GetPartialVolumeAndArea method,
        ' using the triangle enumerator we have established:
        pTinSurface.GetPartialVolumeAndArea nRef, esriPlaneReferenceBelow, pEnumTriangles, , nAreaTemp
        
        ' add this value to the area for the entire feature,
        ' tracking this inside a loop to handle multi-part features:
        nArea = nArea + nAreaTemp
        
      Next
    Else
      ' bUseFeature=false
      nArea = -1
    End If
    
    ' update the field of this feature designated, with the area value:
    Set pRow = pFeature
    pRow.Value(nUpdateField) = nArea
    pRow.Store

    ' stop editing and do not save edits to return the TIN to it's original state:
    pTinEdit.StopEditing False
    
    ' increment to the next feature:
    Set pFeature = pFC.NextFeature
    nFeatures = nFeatures + 1
  Loop

  MsgBox nFeatures & " features were updated with a 3D surface area attribute.", vbInformation, "UpdatePolysWith3DSurfaceArea"

  Exit Sub
  
UpdatePolysWith3DSurfaceArea_ERR:
  On Error Resume Next
  If Not pTinEdit Is Nothing Then pTinEdit.StopEditing False
  MsgBox "Error: " & Err.Description, vbCritical, "UpdatePolysWith3DSurfaceArea"
  
End Sub

'
'  find the layer by name and return its' feature cursor
'
Private Function GetFeatureCursorFromLayer(sLayerName As String) As IFeatureCursor
  On Error GoTo EH
  
  Dim pLayer As ILayer
  Dim pSxDoc As ISxDocument
  Dim pMxDoc As IMxDocument
  Dim pEnumLayers As IEnumLayer
  Dim pFeatClass As IFeatureClass
  Dim pFeatLayer As IFeatureLayer
  Dim pFeatCursor As IFeatureCursor

  Dim pApp As IApplication
  Set pApp = New AppRef
  
  '  get the document
  If TypeOf pApp.Document Is ISxDocument Then
    Set pSxDoc = pApp.Document
    Set pEnumLayers = pSxDoc.Scene.Layers
      
  ElseIf TypeOf pApp.Document Is IMxDocument Then
    Set pMxDoc = pApp.Document
    Set pEnumLayers = pMxDoc.FocusMap.Layers
  End If
      
  If pEnumLayers Is Nothing Then Exit Function
  
  ' find the requested layer:
  Set pLayer = pEnumLayers.Next
  Do While Not pLayer Is Nothing
    If UCase(pLayer.Name) = UCase(sLayerName) Then Exit Do
    Set pLayer = pEnumLayers.Next
  Loop
    
  If pLayer Is Nothing Then
    ' layer not found:
    Exit Function
  End If
    

  ' get the feature cursor:
  Set pFeatLayer = pLayer
  Set pFeatClass = pFeatLayer.FeatureClass
  Set pFeatCursor = pFeatClass.Search(Nothing, False)
    
  ' return:
  Set GetFeatureCursorFromLayer = pFeatCursor
    
  Exit Function
    
EH:

End Function

'
'  accept a layername or index and return the corresponding ILayer
'
Private Function GetLayer(sLayer) As ILayer
  Dim pSxDoc As ISxDocument
  Dim pMxDoc As IMxDocument
  Dim pTOCs As ISxContentsView
  Dim pTOC  As IContentsView
  Dim i As Integer
  Dim pLayers As IEnumLayer
  Dim pLayer As ILayer
  
  On Error GoTo GetLayer_Err

  Dim pApp As IApplication
  Set pApp = New AppRef
  
  If IsNumeric(sLayer) Then
    '  if numeric index:
    If TypeOf pApp.Document Is ISxDocument Then
      Set pSxDoc = pApp.Document
      Set GetLayer = pSxDoc.Scene.Layer(sLayer)
    ElseIf TypeOf pApp.Document Is IMxDocument Then
      Set pMxDoc = pApp.Document
      Set GetLayer = pMxDoc.FocusMap.Layer(sLayer)
      Exit Function
    End If
  
  Else
    '  iterate through document layers looking for a name match:
    If TypeOf pApp.Document Is ISxDocument Then
      Set pSxDoc = pApp.Document
      Set pLayers = pSxDoc.Scene.Layers

      Set pLayer = pLayers.Next
      Do While Not pLayer Is Nothing
        If UCase(sLayer) = UCase(pLayer.Name) Then
          Set GetLayer = pLayer
          Exit Function
        End If
        Set pLayer = pLayers.Next
      Loop
            
    ElseIf TypeOf pApp.Document Is IMxDocument Then
      Set pMxDoc = pApp.Document
      Set pLayers = pMxDoc.FocusMap.Layers

      Set pLayer = pLayers.Next
      Do While Not pLayer Is Nothing
        If UCase(sLayer) = UCase(pLayer.Name) Then
          Set GetLayer = pLayer
          Exit Function
        End If
        Set pLayer = pLayers.Next
      Loop
    End If
  End If
  
  Exit Function
    
GetLayer_Err:
    
End Function


喜欢0 评分0
NinJa
  • 注册日期
  • 发帖数
  • QQ
  • 铜币
  • 威望
  • 贡献值
  • 银元
1楼#
发布于:2003-08-26 10:32
为了先凑够帖子,大家别怪我:)

顶顶顶顶顶顶顶顶顶顶顶顶顶顶顶顶顶顶顶顶顶顶顶顶顶顶顶顶
 顶顶顶顶顶顶 顶顶顶顶顶顶顶顶顶顶 顶顶顶  顶顶顶顶顶
 顶顶顶顶顶顶    顶顶顶顶顶顶顶 顶顶顶顶 顶顶顶顶顶
 顶顶顶顶顶   顶顶顶顶顶顶顶顶顶 顶顶顶顶  顶顶顶顶
 顶顶顶顶  顶顶顶顶顶顶顶顶顶顶顶顶 顶   顶顶 顶顶
 顶顶  顶 顶顶顶顶顶顶顶顶顶顶  顶顶顶 顶顶   顶
 顶顶顶顶顶 顶顶顶顶顶顶顶顶顶   顶顶顶   顶顶顶顶
 顶顶顶顶   顶顶顶顶顶顶顶顶顶顶 顶  顶 顶顶顶顶顶
 顶顶顶顶顶顶顶 顶顶顶顶顶顶顶顶顶  顶  顶   顶顶
 顶顶   顶  顶顶顶顶顶顶顶顶顶 顶顶顶顶  顶顶顶顶
 顶 顶顶   顶顶顶顶顶顶顶顶   顶   顶 顶顶顶顶
 顶顶顶顶顶  顶顶顶顶顶顶顶顶 顶  顶顶 顶 顶顶顶顶
 顶顶顶顶    顶顶顶顶顶顶 顶顶 顶顶顶顶顶 顶顶顶顶
 顶顶顶  顶顶   顶顶顶顶 顶  顶顶 顶顶 顶顶顶顶
 顶   顶顶顶顶    顶顶顶顶顶 顶顶顶   顶顶顶顶
 顶顶顶顶顶顶顶顶顶顶顶顶顶顶顶顶顶顶顶顶顶顶顶顶顶顶顶顶
举报 回复(0) 喜欢(0)     评分
sulh-000
路人甲
路人甲
  • 注册日期2005-03-10
  • 发帖数70
  • QQ
  • 铜币242枚
  • 威望0点
  • 贡献值0点
  • 银元0个
2楼#
发布于:2005-04-04 03:50
<img src="images/post/smile/dvbbs/em01.gif" />
举报 回复(0) 喜欢(0)     评分
study/gis8298
路人甲
路人甲
  • 注册日期2004-09-18
  • 发帖数20
  • QQ
  • 铜币65枚
  • 威望0点
  • 贡献值0点
  • 银元0个
3楼#
发布于:2005-04-17 19:28
<img src="images/post/smile/dvbbs/em02.gif" />
举报 回复(0) 喜欢(0)     评分
eagling
路人甲
路人甲
  • 注册日期2004-10-21
  • 发帖数126
  • QQ
  • 铜币540枚
  • 威望0点
  • 贡献值0点
  • 银元0个
4楼#
发布于:2005-04-20 14:26
<P>看不动啊</P>
中国GIS的未来,要靠GIS农民创造.
举报 回复(0) 喜欢(0)     评分
zhousky
论坛版主
论坛版主
  • 注册日期2003-08-01
  • 发帖数281
  • QQ
  • 铜币1027枚
  • 威望3点
  • 贡献值0点
  • 银元0个
5楼#
发布于:2005-04-20 15:14
代码太长了
不要看我噢
举报 回复(0) 喜欢(0)     评分
nxy_918
路人甲
路人甲
  • 注册日期2003-09-15
  • 发帖数74
  • QQ
  • 铜币325枚
  • 威望0点
  • 贡献值0点
  • 银元0个
6楼#
发布于:2005-04-22 08:40
<P>发代码,总得说明一下,把人家的注释给汉化一下也好,在怎么地要说一下功能,这样贴一堆代码,有点不负责任哦</P>
举报 回复(0) 喜欢(0)     评分
ynkm
路人甲
路人甲
  • 注册日期2004-05-26
  • 发帖数264
  • QQ
  • 铜币27枚
  • 威望0点
  • 贡献值0点
  • 银元0个
7楼#
发布于:2005-04-26 13:38
???
举报 回复(0) 喜欢(0)     评分
goodgoodstudy
路人甲
路人甲
  • 注册日期2005-04-29
  • 发帖数13
  • QQ
  • 铜币140枚
  • 威望0点
  • 贡献值0点
  • 银元0个
8楼#
发布于:2005-04-29 16:50
我晕<img src="images/post/smile/dvbbs/em01.gif" /><img src="images/post/smile/dvbbs/em02.gif" />
举报 回复(0) 喜欢(0)     评分
游客

返回顶部