[原创]自己做插值程序，想和大家分享

楼主^#

更多发布于：2005-09-02 20:33

'''''by kisssy
'''''Kriging 以克里金插值为例
''''strName1:A string that represents your FeatureClass Path
''''strName2:A string that represents the FeatureClass Name
''''sFieldName:the field for Interpolation Public Function Kriging(strName1 As String, strName2 As String,sFieldName As String) As IRasterLayer '克里金 '''''''''''''''''''''''''''''''''操作符 Dim pInterpolationOp As IInterpolationOp Set pInterpolationOp = New RasterInterpolationOp '''''''''''''''''''''''''''''''''操作环境 Dim pEnv As IRasterAnalysisEnvironment Set pEnv = pInterpolationOp
 '''''add shape for setting mask , this is optional Dim pFlayer As IFeatureLayer Set pFlayer = addShp(App.Path + "\Data\BaseDB\JiChuDiLi\city", "大同") 'your shp 'path Dim pGeoDB As IGeoDataset Set pGeoDB = pFlayer.FeatureClass
 Dim pEnvelop As IEnvelope Set pEnvelop = pGeoDB.Extent pEnv.SetExtent esriRasterEnvValue, pEnvelop Set pEnv.Mask = pGeoDB
 ''''set cell size pEnv.SetCellSize esriRasterEnvValue, 600 '600:cellsize
 '''''''''''''''''''''''''''''''''''设置FeatureClassDescriptor Dim pFClass As IFeatureClass Set pFClass = OpenFC2(strName1, strName2) Dim pFDescriptor As IFeatureClassDescriptor Set pFDescriptor = New FeatureClassDescriptor pFDescriptor.Create pFClass, Nothing, sFieldName '''''''''''''''''''''''''''''''''设置搜索半径 Dim pRadius As IRasterRadius Set pRadius = New RasterRadius pRadius.SetVariable 12 'your variant ''''''''''''''''''''''''''''''调用Kriging Dim pOutputRaster As IRaster Set pOutputRaster = pInterpolationOp.Krige(pFDescriptor, esriGeoAnalysisSphericalSemiVariogram, pRadius, True) 
'esriGeoAnalysisSphericalSemiVariogram is 'esriGeoAnalysisSemiVariogramEnum
 '''''''''''''''''''''''''''''输出Raster
 Dim pOutRasLayer As IRasterLayer Set pOutRasLayer = New RasterLayer pOutRasLayer.CreateFromRaster pOutputRaster '''''''''''''''''''''''''着色 UsingRasterClassifyColorRampRenderer pOutRasLayer Set Kriging = pOutRasLayer
End Function
Public Function addShp(strPath As String, strFcname As String) As IFeatureLayer ''''Open WorkSpace Dim myFWKS As IFeatureWorkspace Dim myWKSF As IWorkspaceFactory Set myWKSF = New ShapefileWorkspaceFactory Set myFWKS = myWKSF.OpenFromFile(strPath, 0) If Not myFWKS Is Nothing Then ''''Open Dim myFC As IFeatureClass Set myFC = myFWKS.OpenFeatureClass(strFcname) Dim myDS As IDataset Set myDS = myFC Dim myFLayer As IFeatureLayer Set myFLayer = New FeatureLayer Set myFLayer.FeatureClass = myFC myFLayer.Name = myDS.Name Set addShp = myFLayer End If End Function
Public Function OpenFC2(strPath As String, strFcname As String) As IFeatureClass ''''Open WorkSpace Dim myFWKS As IFeatureWorkspace Dim myWKSF As IWorkspaceFactory Set myWKSF = New ShapefileWorkspaceFactory Set myFWKS = myWKSF.OpenFromFile(strPath, 0) If Not myFWKS Is Nothing Then ''''Open
 Set OpenFC2 = myFWKS.OpenFeatureClass(strFcname)
 End If End Function
Public Sub UsingRasterClassifyColorRampRenderer(pRlayer As IRasterLayer)
 ' ' We're going to create StatsHistogram Dim pRaster As IRaster Set pRaster = pRlayer.Raster Dim pStatsHist As IStatsHistogram Set pStatsHist = New StatsHistogram Dim pCalStatsHist As IRasterCalcStatsHistogram Set pCalStatsHist = New RasterCalcStatsHistogram pCalStatsHist.ComputeFromRaster pRaster, 0, pStatsHist
 ' ' and then classify this data
 Dim pClassify As IClassify Set pClassify = New EqualInterval Dim pClassMaxMin As IClassifyMinMax Set pClassMaxMin = pClassify pClassMaxMin.Maximum = pStatsHist.Max pClassMaxMin.Minimum = pStatsHist.Min Dim Classes() As Double Dim ClassesCount As Long Dim numDesiredClasses As Long 'pClassify.Classify numDesiredClasses pClassify.Classify 8 'class count Classes = pClassify.ClassBreaks ClassesCount = UBound(Classes)
 'Create classfy renderer and QI RasterRenderer interface Dim pClassRen As IRasterClassifyColorRampRenderer Set pClassRen = New RasterClassifyColorRampRenderer Dim pRasRen As IRasterRenderer Set pRasRen = pClassRen
 'Set raster for the render and update Set pRasRen.Raster = pRaster pClassRen.ClassCount = ClassesCount pRasRen.Update
 'Create a color ramp to use Dim pRamp As IAlgorithmicColorRamp Dim pColor As IColor Set pColor = New RgbColor Set pRamp = New AlgorithmicColorRamp pRamp.Size = ClassesCount
 pColor.RGB = RGB(0, 0, 255) 'your color pRamp.FromColor = pColor pColor.RGB = RGB(255, 0, 0) pRamp.ToColor = pColor pRamp.Algorithm = 1 pRamp.CreateRamp True 'Create symbol for the classes Dim pFSymbol As IFillSymbol Set pFSymbol = New SimpleFillSymbol
 'loop through the classes and apply the color and label Dim i As Integer For i = 0 To pClassRen.ClassCount - 1 pFSymbol.Color = pRamp.Color(i) pClassRen.Symbol(i) = pFSymbol pClassRen.Break(i) = Classes(i + 1) Next i
 'Update the renderer and plug into layer pRasRen.Update Set pRlayer.Renderer = pClassRen
 Set pRaster = Nothing Set pRasRen = Nothing Set pClassRen = Nothing Set pRamp = Nothing Set pFSymbol = Nothing End Sub

喜欢0 评分0

个人专栏: https://zhuanlan.zhihu.com/c_165676639

举报回复

xc_zhao 路人甲注册日期2007-11-22 发帖数5 QQ 铜币115枚威望0点贡献值0点银元0个加关注写私信	1楼^# 发布于：2008-12-17 11:49 <P><img src="images/post/smile/dvbbs/em02.gif" /> 厉害</P>
	举报回复(0) 喜欢(0) 评分

coverage 路人甲注册日期2008-07-03 发帖数14 QQ 铜币180枚威望0点贡献值0点银元0个加关注写私信	2楼^# 发布于：2008-09-04 14:10 <P>非常有用！</P><img src="images/post/smile/dvbbs/em01.gif" /><img src="images/post/smile/dvbbs/em02.gif" />
	举报回复(0) 喜欢(0) 评分

license

路人甲

注册日期2003-08-20
发帖数235
QQ33281522
铜币366枚
威望0点
贡献值0点
银元0个

加关注写私信

3楼^#

发布于：2006-12-01 00:44

ding

Gis的小石块 QICQ:33281522 EMAIL:license@vip.sina.com GIS的麦田守望者，希望和大家交流。〓〓〓〓〓〓〓〓〓〓 GISEMPIRE 〓〓灌水★波菜〓〓专用章〓〓〓〓〓〓〓〓〓〓

举报回复喜欢评分

jone1017 路人甲注册日期2005-05-05 发帖数5 QQ 铜币129枚威望0点贡献值0点银元0个加关注写私信	4楼^# 发布于：2006-11-27 19:39 <P>顶！！！！！111</P>
	举报回复(0) 喜欢(0) 评分

liaokobe 路人甲注册日期2006-01-17 发帖数36 QQ 铜币52枚威望0点贡献值0点银元0个加关注写私信	5楼^# 发布于：2006-03-25 11:03 cellsize的值是怎么计算出来的？就好像ArcMap的插值操作里面载入一个点集，自动就有一个cellsize数据，程序中如何获取这个数据呢？
	举报回复(0) 喜欢(0) 评分

Yoyozwf

路人甲

注册日期2006-02-15
发帖数39
QQ
铜币207枚
威望0点
贡献值0点
银元0个

加关注写私信

6楼^#

发布于：2006-02-27 08:35

我正在研究vc的算法，不过我用的是IDW方法生成的Raster，由于用的是对话框设置的一些属性，所以想整理下再贴上来，先贴一段检验license的代码，可以直接用的
BOOL GetSpatialAnalystLicense() { HRESULT hr; IUIDPtr ipUID(CLSID_UID); hr = ipUID->put_Value(CComVariant("esricore.SAExtension.1")); if(FAILED(hr)) return false;
 IExtensionManagerAdminPtr ipExtManAdmResult(CLSID_ExtensionManager);
 hr = ipExtManAdmResult->AddExtension(ipUID,0); if(FAILED(hr)) return false;
 IExtensionManagerPtr ipExtManSpa(ipExtManAdmResult); if(ipExtManSpa == 0) return false;
 IExtensionPtr ipExt; hr = ipExtManSpa->FindExtension(CComVariant("Spatial Analyst"),;ipExt); if(FAILED(hr)) return false;
 IExtensionConfigPtr ipExtConfig(ipExt); if(ipExtConfig == 0) return false;
 esriExtensionState esriSpaExt; hr = ipExtConfig->get_State(;esriSpaExt); if(FAILED(hr)) return false;
 if(esriSpaExt!= esriESUnavailable) { hr = ipExtConfig->put_State(esriESEnabled); if(FAILED(hr)) return false; else return true;
 }
 else { MessageBox("no license!"); return false; }
}

举报回复喜欢评分

Yoyozwf 路人甲注册日期2006-02-15 发帖数39 QQ 铜币207枚威望0点贡献值0点银元0个加关注写私信	7楼^# 发布于：2006-02-16 23:22 <P>支持！</P>
	举报回复(0) 喜欢(0) 评分

wavvylia 路人甲注册日期2003-07-28 发帖数384 QQ 铜币555枚威望0点贡献值0点银元0个加关注写私信	8楼^# 发布于：2005-11-03 15:26 <P>哪种语言不重要，关键是算法和思路。</P>
	举报回复(0) 喜欢(0) 评分

sulin 路人甲注册日期2003-07-28 发帖数158 QQ 铜币501枚威望0点贡献值0点银元0个加关注写私信	9楼^# 发布于：2005-10-17 10:08 <P>大家都是VB代码的多，又没有VC++代码的，希望大侠贴出来，参考参考！</P> <P>好贴，狂顶！</P><img src="images/post/smile/dvbbs/em05.gif" />
	举报回复(0) 喜欢(0) 评分

您需要登录后才可以回帖，登录或者注册

[原创]自己做插值程序，想和大家分享

最新喜欢：