;==========================================================
pro buildshp, maskimage, outlineshpfile, needrange
;trick
maskimage = congrid(maskimage,800,800,/interp,/center)
window, /pixmap, xsize=800, ysize=800
contour, maskimage, path_xy=xy, path_info=info, $
nlevels=1, closed=1, $
xmargin=[0,0], ymargin=[0,0], xstyle=4, ystyle=4
isoShp = obj_new('IDLffShape', outlineshpfile, /update, entity_type=5)
;
; 添加属性
;
isoShp->IDLffShape::AddAttribute, 'LEVEL', 3, 20
value = info.value
uniqV = value[uniq(value, sort(value))]
for i=0, n_elements(UniqV)-1 do begin
Index = where(info.VALUE eq UniqV[i])
x = 0.0
y = 0.0
n_parts = n_elements(Index)
parts = 0
n_vertices = 0
for j=0, n_parts-1 do begin
Pos = index[j]
LL = double(xy[*, info[Pos].offset:(info[Pos].offset+info[Pos].N-1)])
LL[0,*] = LL[0,*] * (needrange[2]-needrange[0]) + needrange[0]
LL[1,*] = LL[1,*] * (needrange[3]-needrange[1]) + needrange[1]
x = [x, transpose(LL[0,*]), (transpose(LL[0,*]))[0]]
y = [y, transpose(LL[1,*]), (transpose(LL[1,*]))[0]]
parts = [parts, parts[j]+n_elements(LL[0,*])+1]
n_vertices = n_vertices + n_elements(LL[0,*]) + 1
endfor
vertices = dblarr(2, n_vertices)
vertices[0,*] = x[1:*]
vertices[1,*] = y[1:*]
parts = parts[0:(n_elements(parts)-2)]
;
; Create structure for new entity.
;
entNew = {IDL_SHAPE_ENTITY}
;
; Define the values for the new entity.
;
entNew.SHAPE_TYPE = 5L
entNew.N_VERTICES = long(n_vertices)
entNew.VERTICES = ptr_new(vertices)
entNew.N_PARTS = long(n_elements(parts))
entNew.PARTS = ptr_new(parts)
isoShp->IDLffShape::PutEntity, entNew
attrNew = IsoShp->IDLffShape::GetAttributes(/ATTRIBUTE_STRUCTURE)
attrNew.ATTRIBUTE_0 = long(i)
isoShp->IDLffShape::SetAttributes, i, attrNew
isoShp->DestroyEntity, entNew
endfor
isoShp->Close
obj_destroy, isoShp
end
;==========================================================
function buildmask, shpfile, needrange
;xsize, ysize可用来控制输出的精度,尺寸上限受系统的限制
xsize = (needrange[2]-needrange[0])*100 < 2000
ysize = (needrange[3]-needrange[1])*100 < 2000
; 创建背板，图形在背板中绘制
window, 1, /pixmap, xsize=xsize, ysize=ysize
wset, 1
oshape = obj_new('IDLffShape', shpfile)
oshape->getproperty, n_entities=nentity
for i=0, nentity-1 do begin
ent = oshape->IDLffShape::getentity(i)
if ptr_valid(ent.parts) then begin
cuts = [*ent.parts, ent.n_vertices]
for j=0, ent.n_parts-1 do begin
x = (*ent.vertices)[0,cuts[j]:cuts[j+1]-1]
y = (*ent.vertices)[1,cuts[j]:cuts[j+1]-1]
x = (x-needrange[0])/(needrange[2]-needrange[0])*xsize
y = (y-needrange[1])/(needrange[3]-needrange[1])*ysize
polyfill, x, y, /device
endfor
endif
oshape->destroyentity, ent
endfor
obj_destroy, oshape
; 拷贝背板中的图像
maskimage = tvrd()
return, maskimage
end
;==========================================================
pro shpoutline
;
;=====STEP 1=====
;在pixmap中绘制所需区域的填充图
;
;源shapefile文件
shpfile = 'C:\RSI\IDL63\resource\maps\shape\cntry02.shp'
;所需区域的范围
needrange = double([105.,0.,122.,23.])
maskimage = buildmask(shpfile, needrange)
;
;=====STEP 2=====
;使用IDL提供的contour方法绘制填充图的等值线，再将等值线保存为shapefile格式。
;
;输出的shapefile文件
outlineshpfile = 'd:\outline.shp'
buildshp, maskimage, outlineshpfile, needrange
end
