最新AE二次开发进行DEM表面积求算.docx
《最新AE二次开发进行DEM表面积求算.docx》由会员分享，可在线阅读，更多相关《最新AE二次开发进行DEM表面积求算.docx（11页珍藏版）》请在冰豆网上搜索。
最新AE二次开发进行DEM表面积求算
AE二次开发进行DEM表面积求算
四、算法代码
此后为自定义类中的代码
（1）栅格数据提取：
采用二维数组对栅格数据进行结果存储，代码为arcengine帮助文档内部搜索出然后修正的
publicdouble[,]ReadWriteRawBlocks（IRasterDataset2rasterDs）
{
//Createaraster.
IRaster2raster2=rasterDs.CreateFullRaster（）asIRaster2;
//Createarastercursorwithasystem-optimizedpixelblocksizebypassinganull.
IRasterCursorrasterCursor=raster2.CreateCursorEx（null）;
//UsetheIRasterEditinterface.
IRasterEditrasterEdit=raster2asIRasterEdit;
//Loopthrougheachbandandpixelblock.
IRasterBandCollectionbands=rasterDsasIRasterBandCollection;
IPixelBlock3pixelblock3=null;
IRawBlocks rawBlocks=（IRawBlocks）bands.Item（0）;
IRasterInfo rasInfo=rawBlocks.RasterInfo;
longblockwidth=0;
longblockheight=0;
System.Arraypixels;
IPnttlc=null;
objectv;
double[]v1=newdouble[1000000];
intm=0;
intlenght=0;
double[,]value=newdouble[rasInfo.Width,rasInfo.Height];
longbandCount=bands.Count;
do
{
pixelblock3=rasterCursor.PixelBlockasIPixelBlock3;
blockwidth=pixelblock3.Width;
blockheight=pixelblock3.Height;
pixelblock3.Mask（255）;
//value=newdouble[blockwidth,blockheight];
for（intk=0;k {
//Getthepixelarray.
pixels=（System.Array）pixelblock3.get_PixelData（k）;
for（longi=0;i {
for（longj=0;j {
//Getthepixelvalue
v=pixels.GetValue（i,j）;
//Dosomethingwiththevalue.
//此处的数据提取是以64*64单元提取，需要主意顺序，该段代码前一讲中已说明。
v1[m]=Convert.ToDouble（v）;
lenght=m;
m++;
}
}
//Setthepixelarraytothepixelblock.
pixelblock3.set_PixelData（k,pixels）;
}
//Writebacktotheraster.
tlc=rasterCursor.TopLeft;
rasterEdit.Write（tlc,（IPixelBlock）pixelblock3）;
}
while（rasterCursor.Next（）==true）;
System.Runtime.InteropServices.Marshal.ReleaseComObject（rasterEdit）;
m=0;
//此处的循环为了将提取的数据进行重新排序
while（m {
for（inti=0;i for（intj=0;j for（intb=0;b<64;b++）
for（intd=0;d<64;d++）
if（i+d {
value[j+b,i+d]=v1[m];
m++;
}
}
returnvalue;
}
（2）三角型面积计算，cellx，celly为栅格每一格所代表的实际长度
publicdoubleAreaBase（doublecellX,doublecellY,doublecellZ1,doublecellZ2,doublecellZ3）
{
doublevalue1,value2,value3;
value1=Math.Sqrt（cellX*cellX+（cellZ1-cellZ2）*（cellZ1-cellZ2））;
value2=Math.Sqrt（cellX*cellX+（cellZ3-cellZ2）*（cellZ3-cellZ2））;
value3=Math.Sqrt（cellX*cellX+cellY*cellY+（cellZ3-cellZ1）*（cellZ3-cellZ1））;
doublelength=value1+value2+value3;
return Math.Sqrt（length*（length- value1）*（length- value2）*（length- value3））;
}
（3）总体面积计算：
此处采用一个sbyte三维数组，将每个栅格格子分成八个小三角形，目的是为了记录该快小三角形是否被之前相邻的格子使用过。
以免造成面积重复计算
publicstring AreaCalculate（）
{
try
{
IRasterBandCollectionrbcollection=rasterasIRasterBandCollection;
IRasterDatasetrasterDataset1=newRasterDatasetClass（）;
rasterDataset1=rbcollection.Item（0）.RasterDatasetasIRasterDataset;
IRasterDataset2rasterDataset2=rasterDataset1asIRasterDataset2;
doublecellX,cellY;
intXnum,Ynum;
RasterInforasterInfo=newRasterInfo（）;
rasterInfo.GetRasterInfo（raster）;
Xnum=rasterInfo.Height;
Ynum=rasterInfo.Width;
cellX=（rasterInfo.xmax-rasterInfo.xmin）/Xnum;
cellY=（rasterInfo.ymax-rasterInfo.ymin）/Ynum;
double[,]value=newdouble[Ynum,Xnum];
sbyte[,,]bool1=newsbyte[Ynum,Xnum,8];
for（inti=0;i for（intj=0;j for（intk=0;k<8;k++）
bool1[i,j,k]=0;
value=ReadWriteRawBlocks（rasterDataset2）;
doublearea=0;
doublesarea=0;
//returnvalue[0,0];
for（inti=1;i {
for（intj=1;j {
if（value[i,j]!
=-32768）
{
if（value[i,j-1]!
=-32768&&value[i+1,j]!
=-32768）
{
if（bool1[i,j,4]==0&&bool1[i,j,5]==0&&bool1[i,j-1,7]==0&&bool1[i+1,j,2]==0）
{
sarea=AreaBase（cellX,cellY,value[i,j-1],value[i,j],value[i+1,j]）;
area=area+sarea;
bool1[i,j,4]=1;
bool1[i,j,5]=1;
bool1[i,j-1,7]=1;
bool1[i+1,j,2]=1;
}
}
if（value[i+1,j]!
=-32768&&value[i+1,j-1]!
=-32768）
{
if（bool1[i,j,5]==0&&bool1[i+1,j,2]==0&&bool1[i+1,j,3]==0&&bool1[i+1,j-1,0]==0）
{
sarea=AreaBase（cellX,cellY,value[i,j],value[i+1,j],value[i+1,j-1]）;
area=area+sarea;
bool1[i,j,5]=1;
bool1[i+1,j,2]=1;
bool1[i+1,j,3]=1;
bool1[i+1,j-1,0]=1;
}
}
if（value[i+1,j]!
=-32768&&value[i+1,j+1]!
=-32768）
{
if（bool1[i,j,6]==0&&bool1[i+1,j,0]==0&&bool1[i+1,j,1]==0&&bool1[i+1,j+1,3]==0）
{
sarea=AreaBase（cellX,cellY,value[i,j],value[i+1,j],value[i+1,j+1]）;
area=area+sarea;
bool1[i,j,6]=1;
bool1[i+1,j,0]=1;
bool1[i+1,j,1]=1;
bool1[i+1,j+1,3]=1;
}
}
if（value[i,j+1]!
=-32768&&value[i+1,j]!
=-32768）
{
if（bool1[i,j,6]==0&&bool1[i,j,7]==0&&bool1[i,j+1,4]==0&&bool1[i+1,j,1]==0）
{
sarea=AreaBase（cellX,cellY,value[i,j+1],value[i,j],value[i+1,j]）;
area=area+sarea;
bool1[i,j,6]=1;
bool1[i,j,7]=1;
bool1[i,j+1,4]=1;
bool1[i+1,j,1]=1;
}
}
}
}
}
doublemarea=cellX*cellY/8;
for（inti=0;i for（intj=0;j {
if（value[i,j]!
=-32768）
{
for（intk=0;k<8;k++）
{
if（bool1[i,j,k]==0）
area=area+marea;
}
}
}
returnarea.ToString（）;
}
catch（Exceptionex）
{
returnex.ToString（）;
}
}
以下为主函数代码，写在mapcontrol的mousedown事件中
（4）栅格裁剪：
首先交互式绘制矢量多面性，采用IExtractionOp接口对原始栅格进行裁剪提取。
IGeometrygeometry=null;
geometry=axMapControl1.TrackPolygon（）;
drawMapShape1（geometry,200,50,0）;
IPolygonpPolygon=geometryasIPolygon;
IExtractionOpextraction=newRasterExtractionOpClass（）;
升级会员 