科学网

 找回密码
  注册

tag 标签: 等值线

相关帖子

版块 作者 回复/查看 最后发表

没有相关内容

相关日志

等离子体粒子模拟及其应用课程大作业暨2020年第4次组会报告(附教学视频)
YUNJU 2020-5-26 00:07
目录 程序变量说明 nss:两种成分的电子 np:粒子数目 idimv:每个电子有两个坐标,三个速度 nx,ny:分别是x,y方向的网格点 nloop:计算的总步数 dt:时间步长 part(idimv,np,nss):所有电子坐标和速度的数组 dns(nss):两种电子成分的权重 part(idimv,np):所有离子坐标和速度的数组 vbe(idimp),vte(idimp):电子初始的平均速度和热速度 vbi(idimp),vti(idimp):离子初始的平均速度和热速度 fx(nxv,ny),fy(nxv,ny):网格上电场的x,y分量(实空间) fxc(nxv,ny),fyc(nxv,ny):网格上电场的x,y分量(复空间) qs(nxv,ny,nss):两种电子成分在网格上的电荷 q(nxv,ny):网格上的电荷密度(实空间) qc(nxv,ny):网格上的电荷密度(复空间) we:通过电势算的电场能量 eke:电子动能 eki:粒子动能 wt:系统中的总能量 fx:电场在x方向的分量 fy:电场在y方向的分量 wx:x方向上的电场能量 wy:y方向上的电场能量 验证数值不稳定性 网格自加热效应验证 束流不稳定性结果(附MATLAB代码) clc;clear all;load D:\\Fortran\\shuliu\\ee.dat; loop=15;Cxy=256;qq=zeros(Cxy,Cxy,loop); ex=zeros(Cxy,Cxy,loop);ey=zeros(Cxy,Cxy,loop); for i=1:loop for j=1:Cxy for k=1:Cxy temp=Cxy*(j-1)+k+Cxy*Cxy*(i-1);qq(j,k,i)=ee(temp,3); ex(j,k,i)=ee(temp,4);ey(j,k,i)=ee(temp,5); end end end x=1:Cxy;y=1:Cxy; =meshgrid(x,y) figure(1) for i=1:loop subplot(5,3,i);contour(X,Y,qq(:,:,i)) end figure(2) for i=1:loop subplot(5,3,i);contour(X,Y,ex(:,:,i)) end figure(3) for i=1:loop subplot(5,3,i);contour(X,Y,ey(:,:,i)) end 附教学视频(有两个大于4G的视频未上传): https://pan.baidu.com/s/1caKhAKCmdKRibdgF2aJJ0A
个人分类: 科研笔记|121 次阅读|1 个评论
全球2020年春季新冠病毒感染扩散势态等值线图剖析
bqzhu 2020-4-27 19:55
今年春季将过,夏季将来,全球新冠病毒的扩散引起了人们对生命安全与未来经济发展的担忧与恐慌 。 根据相关全球公布的新冠病毒确诊感染人数数据库和人口数据库资料,进行了 全球 2020 年春季新冠病毒感染扩散势态等值线 数据处理与作图,数据为确诊感染人数据占人口比例的千分数。 从等值线图剖析 , 可以得到以下几点重要认识 : (1) 大规模 全球新冠病毒感染率強与超強区只出现在北大西洋周边国家与地区,並向东西方向扩散並逐渐减弱。 西欧強与超強感染出现在西 班牙,英国,意大利,卢森堡,比利时,瑞士,爱尔兰,法国等国 ( 确诊 感染率为 2-6.5) 。德国与奥地利减弱为重感染区 (1.5-2) ,东欧地区 一 在 0.5 左右,中东地区为 0.5-1 左右。而到中亚地区则减至 0.2 以下。北欧近北大西洋的瑞典,挪威 确诊 感染率是偏离大西洋的芬兰的 2-3 倍。 美国超強感染区出现在东北角大西洋沿岸,包括 纽 约、新泽西、麻萨诸塞 、 康涅狄克 、 罗得島 等 州 ( 确诊 感染率为 6-15) 。 从西北到东南太平洋沿岸则明显梯度下降 ( 科罗拉多为 2 ,加里福尼亚为 1) ,向西变化也是如此 ( 密歇根在 3.5 左右,华盛顿州在 1.5 左右 ) 。美国大西洋岸由北向南疫情变化则起伏不定,南部路易斯安那州又成为超強感染区 ( 感染率 7) 。 加拿大疫情的強感染区出现在大西洋边的魁北克省 ( 确诊感染率 2.5) ,占全国确诊感染人数一半以上。而到太平洋边的英属哥伦比亚省确诊感染率降到 0. 5 以下。 位于北大西洋北部北纬 65 度左右的冰島是一个人居和外来人均很少且分散的地方,一个想要碰到个人都难的地方,但恰是一个超強感染区 ( 确诊感染率 5.2) 。位于北纬 31 度的百慕大群島更是一个人居稀少,只有少数冒险家敢去的地方也是重感染区 ( 感染率近以 2) ,是全球小島国之最。 (2) 全球连成片的重感染 - 超強感染区明显受到经纬度制约。北纬 30 度左右是南向的控制线, N35 度以南感染率急剧下降。希腊、意大利南部西西里島、塞浦路斯等地区为弱感染区。非洲北部目前确诊感染率在 十 万分几,可能与检测量不够有关,但即使确诊感染数增加一个数量级,也仍然是弱感染区,而且是水涨船高,欧美会增长更多。因此不会改变图的势态。 N30 度线左右是地球生态环境变化和全球突发 、 罕见事件出现的重要控制线。虽然机理不很清楚,但一直受到古今中外科学家,探險家考古学家的关注,如神秘的百慕大三角位于 N30 度线。因此冬春季节重至超强新冠感染区主要出现在北半球温带地区,向寒带地区过渡则趋于减弱。 欧美重至超強感染区向东扩散至东径 70 度线总体上趋于终止。青藏高原与中亚高原无论是阻止人流传布还是自然传布都起着天然屏障作用。 (3)E70 度以东的中国、东南亚 、 蒙古 、 朝鲜等地总体上为疫情关注区;中国除湖 北 外的各省市自治区确诊感染率低于十万分之 3 (0.03) ,平均为 十 万分之 1.2(0.012) ,东南亚地区除新加坡 、 马来西亚外,也均小于十万分之 2(0.02) 。西藏 、 新疆 、 蒙古 、 朝鲜总体上无疫情。 但在这片关注区出现了两个孤島。一个是新加坡,为外来输入。由于隔离措施迟缓,亮起了红灯,为強感染区 ( 确诊感染率 2) 。另一个是武汉,由于早发现早隔离,使这个孤島变得非常更小,湖北被控制为重感染区 ( 确诊感染率为 1.15) 。武汉地区没有出现疫情扩散面,因此不存在源区性质的自然传布,明显是境外输入型。韩国疫情也具孤島性质,总体上为弱感染区 ( 确诊感染率为 0.2 左右 ) ,为外来输入。 (4) 病毒传布除 了 重视与控制人流传播外,应充分重视源区性质的自然传播的调查与控制。特別在病毒产生的早期阶,自然传布可能更为重要。有以下一些现象提供了自然传播的信息与思考 : (a) 为什么冬春季节新冠病毒在 北 半球溫带地区大规模流行扩散?这是个自然因素,不是人为因素。是否存在受北半球温带季风影响的气溶胶传播?意大利科学家已发现气溶胶中存在病毒。不 少 学者认为病毒与生命起源有关,发育的早期阶段是非生命性质的有机分子团,因此不存生与死,可以被气溶胶携带远距离传播。 (b) 为什么在冰島 、 魁北克、百慕大这些人烟稀少的地方会出现強至超強感染区?显然不是人流传播。 据 deCODE 基因公司 资料 , 在冰岛地区 对确诊患者进行 基因 检测,发现 有 40 种变异病毒 ,并发现首例二次感染者, 体内存在 2 种新冠病毒,其中一种是原始的变体 ,无证状感染者比例高达 4 8% ,显然病毒在冰島已存在长时间的流传。 (c) 为什么北大西洋周边地区病毒感染传播形成了扩散趋势面和梯度变化 ? 自然界物质与能量传布服从扩散方程描述,受物质浓度和能量密度梯度制约。大众都明白自然界水只能从高处往低处流,难道西方 的 一些政客不明白这个道理。而人流传播常是孤島式,湍乱式,难于形成有规律的扩散趋势面。 全球 2020 年春季新冠病毒感染扩散势态 等值线 图 图例 : 疫情等值分区(人口 确诊 感染率千分数) : 0.02 关注区 , 0.02-0.1 控制区 , 0.1-0.5 弱感染区 , 0.5-1 感染区 , 1-2 重感染区 , 2-4 强感染区 , 4 超强感染区 . 。
个人分类: 环境灾害|3647 次阅读|0 个评论
【Matlab】如何填充contour等值线的陆地颜色
JerryYe 2016-8-27 09:53
第一步,先用contourf画出相应的图形,比如下图: 第二步,打开colormap的编辑器“Edit--Colormap” 第三步,将数值为“0”的等值线颜色设置为黑色(按需要修改),其他颜色设置为白色 效果如下: 第四步,在上图基础上,添加等值线contour(lon,lat,topo, ,'k','linewidth',2) 效果如上图.
个人分类: Matlab|17407 次阅读|0 个评论
土壤湿度等值线Contourf
YF2015 2016-5-28 20:30
contourf Filled two-dimensional contour plot Syntax contourf(Z) contourf(Z,n) contourf(Z,v) contourf(X,Y,Z) contourf(X,Y,Z,n) contourf(X,Y,Z,v) = contourf(...) Description A filled contour plot displays isolines calculated from matrix Z and fills the areas between the isolines using constant colors. The color of the filled areas depends on the current figure's colormap. contourf(Z) draws a contour plot of matrix Z, where Z is interpreted as heights with respect to a plane. Z must be at least a 2-by-2 matrix. The number of contour lines and the values of the contour lines are chosen automatically. contourf(Z,n) draws a contour plot of matrix Z with n contour levels. contourf(Z,v) draws a contour plot of matrix Z with contour levels at the values specified in vector v. contourf(X,Y,Z) , contourf(X,Y,Z,n) , and contourf(X,Y,Z,v) produce contour plots of Z using X and Y to determine the x - and y -axis limits. When X and Y are matrices, they must be the same size as Z, in which case they specify a surface as surf does. = contourf(...) returns the contour matrix C as calculated by the function contourc and used by clabel, a vector of handles h to patch graphics objects, and a contour matrix CF for the filled areas. Remarks If X or Y is irregularly spaced, contourf calculates contours using a regularly spaced contour grid, then transforms the data to X or Y. Examples Create a filled contour plot of the peaks function. = contourf(peaks(20),10); colormap autumn See Also clabel , contour , contour3 , contourc , quiver Contour Plots for related functions 土壤湿度等值线
个人分类: 论文写作|3489 次阅读|0 个评论
变量A的等值线分布图与其倒数的分布图的关系是这样的!
热度 1 zhangxw 2015-9-29 11:51
变量 A 的等值线分布图与其倒数的分布图的关系是这样的! 张学文, 2015.9.29 最近我计算了新疆各地的水分循环一次的平均天数(不妨以 A 表示它)的等值线的地理分布的草图 http://blog.sciencenet.cn/home.php?mod=spaceuid=2024do=blogid=922824 ,朋友檀成龙建议我再分析该变量的倒数,一年水分循环的次数( B )的分布图。 是的,各地的水分循环平均每年的循环次数也是一个气候变量。分析这个新疆(或者全国 / 全世界)各地一年的水汽的平均循环次数,确实是他人没有做过的一件事。为什么不做一下,何况数据几乎是现成的。 昨天我就动手分析新疆的这种分布图,并且绘出了等值线。今天再看这张一年的水分循环平均次数图,发现它与前面分析的水分循环一次的天数分布图十分类似。在分布图上吐鲁番的高值中心(循环一次需要 302 天)现在反而变成了低值中心(一年循环的次数最少)。 于是我问自己除了等值线的高低中心相反之外,是否等值线的走向也一致? 我考虑了两分钟。。。 我最后明白了,由于这两个变量互为倒数,原来循环一次需要的时间分布图和一年的循环次数图上的等值线应当全部的平行的,而其等值线的梯度方向是相反的。并且高中心都变成了低中心 / 低中心变成了高中心。 这不仅是我分析的结果是如此,而且所有的变量 A 的等值线与 A 的倒数( 1/A )的等值线都应当是如此(条件是变量大于 0 )。这应当是个自然的 数学结论 。而水分循环一周需要的时间与每年水分循环次数仅是这个通用原则的特例! 这样看来,分析水分循环一种分布图就基本上够了(我昨天分析的图也不必公布了)。 看来我画了一张可以不画的分布图,但是我明白了 互为倒数的变量的等值线的地理分布图的关系是: 1. 高 / 低中心是互换的; 2. 等值线是互相平行的; 3. 等值线的梯度方向是相反的; 4. 等值线的密度是不同的。 原来,变量 A 的等值线分布图与其倒数的分布图的关系是这样的。
个人分类: 空中水科学|2327 次阅读|1 个评论
[转载]如何在后处理时显示等值线+标注
Daniel1985 2012-10-7 10:18
原文地址: http://bbs.caetecc.com/thread-550-1-1.html Contour Plot Option--basic--coutour type 选择line 顺便谁能在coutour line 上显示数值,可以共享一下 网上有一种说法,可以在viewport--create Annotation--text 但在这种方法只能解决表面问题,相当于在屏幕上贴了一行文字。
个人分类: ABAQUS/PYTHON/FORTRAN|3556 次阅读|0 个评论
环境信息经验一:等值线程序
热度 11 Talky 2012-5-19 16:04
环境信息经验一:等值线程序 八十年代末期从事环境影响评价工作时,一个困难是没有绘图软件(程序),地图和等浓度线画起来难度大。 1989 年自行开发完成了等值线程序,得到了很多应用。为了提高计算速度,程序是 FORTRAN 写的,计算应该连接的各点坐标,记录在数据文件中。然后通过 BASIC 程序,调用这个数据文件,结合地图,调用绘图仪输出。之后又改写为 C++ 语言模块,放在环境信息系统中应用。原因是地理信息系统( GIS )非常商业化,而环境信息系统常常主要应用 GIS 中的地图功能,“嵌入式”是经常选择。这时 GIS 往往不含等值线作图部分,这时自己开发的等值线模块就有用了。 这里提供思路和源程序,供有兴趣的青年人应用。 近几年没有再仔细调研过等值线程序的国内外文献,以前做过,感觉我的方法是比较好的。想过发文章,一直懒,最终也没有写好,没有发表过。曾经想称为:“ 全局定位型矩形网格等值线算法" 思想:如图 1 ,如果已经有计算好的网格浓度,用颜色深浅表示。 图中有 4 种不同深浅的颜色(应用计算机很容易分辨的)。如果等值线就按这 4 个值,每个值有几条线,在哪里,已经一目了然了。因此可以利用这一步,可能联结的每条等值线上的点坐标,对应的网格点“点对”。然后按这些记录下来的“点对”数据,一条一条地排序判断,连接次序,计算真正需要的点坐标。 这个算法的优点是,不可能发生连接误差:交叉或遗漏。缺点是要先计算网格点浓度,精度受到网格疏密的限制,因此也受到计算语言软件数组大小的限制。经验是, 81 81 的网格在一般情况下足够了。图中是 23 23 。计算结果如图二: 附源程序: FORTRAN, BASIC, C++ 。 供参考,是 1994 年的,去掉了一些不必要的字符。 DZX-FOR.FOR DZX-BAS.BAS DZX-CPP.CPP
个人分类: 环境信息|7198 次阅读|25 个评论
[转载]C#绘制等值线
热度 1 majian 2012-2-14 16:05
转载自: http://www.4ucode.com/Study/Topic/2104838 http://www.cnblogs.com/think8848/archive/2011/05/04/2036179.html 有一本相关的书:《限定Voronoi网格剖分的理论及应用研究》由北京邮电大学出版社出版。全面介绍了限定Voronoi图的概念、生成技术,采用灵活性更好的带权Delaunay三角/四面体剖分来解决二维/三维限定Voronoi网格剖分的问题,所得到的限定Voronoi网格具有同限定Delaunay三角网格相似的优良性质。建立起了二维/三维限定Voronoi的质量和尺度评价准则,设计了二维/三维限定Voronoi网格的质量和尺度控制的算法。最后,给出了限定Voronoi网格剖分的一些应用实例。 本书可供计算几何、地理信息系统、机器人、通信i石油地质勘探及其相关领域的科研人员及高等学校相关专业师生参考使用。 ======================================== 这两天手头有个项目,需要绘制等值线,本以为是一个很简单的事情,没有想到刚开始就发现竟然无从着手,研究了一个星期,终于把线条画出来了,基本思路是先三角网剖分,然后再等值线追踪,最后绘制;没有对等值线进行光滑处理,示例图中看起来比较光滑是因为取点比较密集,也没有打算进行等值线填色,因为项目中没有这个需求,(而且在我的项目中高程点是网格状分布,而不是离散点,因此我做的三角网剖分简单,但是等值线追踪算法是完全满足离散点要求的)。先上几个效果图: 示例图(黄颜色圆圈代表光源,高程值为光源照度) 图1 图2 图3 等值线标注示意图 效果一:高程值压线了 效果二:高程值在线条下方 效果三:高程值没有按照等值线方向旋转 效果四:在很小的封闭式等值线圈内进行标注 绘制等值线的理论有很多,但是通常大家采用的是三角网剖+等值线追踪,网上有很多论文可以参考: Delaunay三角网剖分算法 Delaunay三角剖分(Delaunay Triangulation)相关知识 三角网格等值线自动生成方法及程序实现 这里 是一个绘制、填充等值线的源代码(《C#实现基于三角网的等值线追踪及填充算法》) 这里 是一些三角网剖分源代码,有各种语言的实现 这里 是一些其他等值线绘制方法 这里 是论文的链接地址,该链接包含很多论文地址。 基本上来说,在《三角网格等值线自动生成方法及程序实现》一文中,提供的理论就足以指导绘制作业了,但是不幸的是,相当多的程序员可能不会一下子就能行成一个比较清晰明确的思路,我就是其中之一,看了好多的文章才明白,噢!原来是这么回事啊。另外还是那句老话:“眼里过千遍,不如手里过一遍”,所谓纸上得来终觉浅,绝知此事要躬行,看算法,看代码,可能比较难以一下子明白其中的原理,很多时候人家用了一个算法(或数据结构),你可能很不以为然,但是当你自已动手时,才发现原来那样用有那样用的妙处(或只能那么用),眼高手低永远是我等程序员的通病。 ok,进入正题 在我的项目中,是计算在一个空间内部,假设排布了n盏灯,那么需要计算出空间内部每一个点上的照度是多少,同时绘制出照度的等值线来,计算照度不在本文计论范围之内,于是我们假设把空间划分成了一个矩阵,矩阵中每个点上的照度值均已获得。于是很自然想到“使用网格绘制等值线”这么一个关键词,发现居然很少有相关资料,于是只能把把关键词修改为“绘制等值线”,于是发现了很多关于绘制等值线方法的论文,但是基本上前面都带了另一个词:三角网。为什么要三角网呢?我们先来看一下使用网格绘制时一种情形 图1 在使用网格绘制等值线时,假定我们是从这个网格的左边起开始查找等值线的走向,则等值线有三个方向可以“出去”,判断和确定等值线的走向将是一个非常困难的事情,而且极有可能出现锯齿状等值线(如下图所示) 图2 而使用三角网为基础来绘制等值线时,则是下图的情况 图3 在这里我有一个问题,到目前为止还没有解决,但是也没有找到解决问题的理论依据,在使用网格时,需要判断等值线的走向,有个论文中提出一个原则:假设等值线的“入口”在左边上(图1A点),则在网格中,如果等值线的“出口”有多个(图1B,C,D点),那么首先要看到等值线到A的趋势更倾向于择哪个点(B,C,D),其次是看哪个点到A点的距离近则取哪个。但是实际上这样做的代价太大了,任何人要去写一段判断趋势的代码,肯定都是个头疼的问题,因此实际上的作法是:在由左向右的追踪过程中,先看哪个点的X值和A点的X值接近则取哪一个,如果有多个点的X值和A点的X值相等,则看哪个点到A点的距离小则取哪一个,如果这个距离也相等,则只能通过趋势判断了(见《 基于规则网格等值线生成算法及其应用 》1.2.2)。实际从程序员的角度来讲,只要你明确表示在程序中会用到某个功能,那么你一年用n+1次和用1次是没有区别的,因此如果你想把程序写的能覆盖所有的可能,判断趋势看起来是不可避免的。 虽然三角形简化了上述这种在网格中追踪等值线的难度,但是并没有完全解决需要判断趋势的问题,为啥咧?使用三角形时只可能有两个“出口”,使用网格中的原则基本上就可以判断取哪个点了,但是理论上在点A所在的边上,仍一个点,到B,和C的X值和距离是相等的,见下图 图4 在图4中,三角形ABC是一个等腰三角形,这就意味着A点到B以及C的距离相等,而且B和C的X值也相等...理论上来讲,如果A是起点,哪还谈什么趋势呢?还要不要让人活了... OK,现在让我们一起默念“阿门”,求助于概率吧,不要出现这种情况,反正我是不打算处理这种情况的,出了Bug再说吧。值的说明的是,在约大多数的论文中,没有提到该如何去处理这情况,有个别论文中指出:“每个三角形上不可能三边都有同值的等值点,另一边上必定有同值的等值点”(见《 如何根据离散点自动绘制等值线(等高线)之 三角形法 》1.b),但是这个说法我持怀疑态度,但我现在不打算处理这种情况(实际上,我在代码中连距离和X值都没有处理)。 有了理论作为依据,接下来的事情就是先把矩阵转换为三角网,怎么转换是个问题:实际上很多时候关于Delaunay的论文都是在假定数据(如灯具照度、高度、降水量等)是离散数据,即不是规则分布的,因此在这些论文中都是在讲如果把这些离散点划分为三角形,当然,能把离散点剖分出三角网,则么规则点就不在话下了。但在我的项目中是不存在这些需求,因此我需要做的就是把矩阵中的每个网格一分为二,切成两个直角三角形,一率从左上角向右下角切,如图所示: 图5 切的有点密,这个可以视具体情况自已定夺。 在这里,需要把三角形、边、点的数据结构给出来: VPoint代表一个点,由X,Y轴坐标及高程值构成 ? public class VPoint { public float X { get ; set ; } public float Y { get ; set ; } public decimal Value { get ; set ; } public override bool Equals( object obj) { if (obj is VPoint) { var tmp = obj as VPoint; return this .X == tmp.X this .Y == tmp.Y; } return false ; } public override int GetHashCode() { return base .GetHashCode(); } } Edge代表一个边,由两个点构成,这里需要稍作解释,在图5中,有一些三角形的一条边是整个三角网的边框,即空间的四个边,开放式的等值的起点和终点必然会落在这些最外层的边上,因此需要对这些最外层边作一个标记,如何标记呢?很显然,三角网中除了这些最外层边,其他所有的边都会同时属于两个三角形,因此,我们的Edge对象有两个属性:Trangle1和Trangle2,这两个属性分别标识这个Edge对象所属的两个三角形,那么这些最外层的Edge对象比较可怜,他们只属于一个Trangle,属性Trangle2的值一定是null,参见Trangle类型。 ? public class Edge { public VPoint P1 { get ; set ; } public VPoint P2 { get ; set ; } public Triangle Triangle1 { get ; set ; } public Triangle Triangle2 { get ; set ; } public bool IsBorder { get { return Triangle2 == null ; } } public override bool Equals( object obj) { if (obj is Edge) { Edge tmp = obj as Edge; return this .P1 == tmp.P1 this .P2 == tmp.P2; } return false ; } public override int GetHashCode() { return base .GetHashCode(); } } Trangle代表一个三角形,由三个边构成,我总是把水平的那条边称为A,竖直的那条边称为B,斜边称为C,其他的属性值以后会用,这里就不一一介绍了, 当然,现在代码还没有进行最终整理,有可能有些属性也许会被干掉,最终只留下必要的。 ? public class Triangle { private Edge a, b, c; /// summary /// 邻边 /// /summary public Edge A { get { return this .a; } set { this .a = value; if ( this .a.Triangle1 == null ) { this .a.Triangle1 = this ; } else { this .a.Triangle2 = this ; } } } /// summary /// 对边 /// /summary public Edge B { get { return this .b; } set { this .b = value; if ( this .b.Triangle1 == null ) { this .b.Triangle1 = this ; } else { this .b.Triangle2 = this ; } } } /// summary /// 斜边 /// /summary public Edge C { get { return this .c; } set { this .c = value; if ( this .c.Triangle1 == null ) { this .c.Triangle1 = this ; } else { this .c.Triangle2 = this ; } } } /// summary /// 是否使用,-1:临时使用; 0:未使用; 1:使用 /// /summary public int UsedStatus { get ; set ; } public Edge BorderEdge { get { if ( this .A.IsBorder) { return this .A; } if ( this .B.IsBorder) { return this .B; } return null ; } } public Edge edges = new Edge ; if (relativeEdge == this .A) { edges = this .B; edges = this .C; } else if (relativeEdge == this .B) { edges = this .A; edges = this .C; } else { edges = this .A; edges = this .B; } return edges; } public VPoint FindPoint( float x, float y) { if ( this .A != null ) { if ( this .A.P1.X == x this .A.P1.Y == y) { return this .A.P1; } if ( this .A.P2.X == x this .A.P2.Y == y) { return this .A.P2; } } if ( this .B != null ) { if ( this .B.P1.X == x this .B.P1.Y == y) { return this .B.P1; } if ( this .B.P2.X == x this .B.P2.Y == y) { return this .B.P2; } } if ( this .C != null ) { if ( this .C.P1.X == x this .C.P1.Y == y) { return this .C.P1; } if ( this .C.P2.X == x this .C.P2.Y == y) { return this .C.P2; } } return null ; } public Edge FindEdge(VPoint p1, VPoint p2) { if ( this .A != null ) { if ( this .A.P1 == p1 this .A.P2 == p2) { return this .A; } } if ( this .B != null ) { if ( this .B.P1 == p1 this .B.P2 == p2) { return this .B; } } if ( this .C != null ) { if ( this .C.P1 == p1 this .C.P2 == p2) { return this .C; } } return null ; } } OK,现在是时候说说我是如何建立一个IListTrangle的吧,即生成三角网列表 一开始,需要解决“有木有”的问题,于是我先建立一个IListTrangle的实例result,然后生成四个VPoint对象,代表一个网格的四个角,然后再生成五个Edge对象,再使用五个Edge生成两个Trangle对象,然后把两个Trangle对象加入到result中,以后每次都去查找result中有没有指定X,Y的VPoint,以及查找result中有没有指定P1和P2的Edge,如果有的话,就分别拿出来去构成新的Edge(对于VPoint来说)或者Trangle(对于Edge来说),如果没有的话就new一个实例出来,然后再通过Trangle对象存到result中,参见Trangle.FindPoint方法和Trangle.FindEdge方法,因为大量使用Linq方法查询数据,造成的结果是45× 27的矩阵,生成三角网需要花费00:00:04.7031250!后来我改进了方法,使用一个二维数组来保存VPoint对象实例(VPoint ),使用一个四维数组来保存Edge对象实例(Edge ),直接用矩阵的序号来获取可能已存在的实例,果然性能得到极大的提升,运行耗时仅00:00:00.0156250。 ? for ( int x = 0; x xCount - step; x += step) { for ( int y = 0; y yCount - step; y += step) { var p1 = matrix ; var p2 = matrix ; var p3 = matrix ; var p4 = matrix ; VPoint tmpP1, tmpP2, tmpP3, tmpP4; Edge edge1, edge2, edge3, edge4, edge5; Triangle triangle1 = new Triangle(), triangle2 = new Triangle(); tmpP1 = this .FindOrCreateNewPoint(tmpMatrix, x, y, p1, zoomFactor); tmpP2 = this .FindOrCreateNewPoint(tmpMatrix, x + 1, y, p2, zoomFactor); tmpP3 = this .FindOrCreateNewPoint(tmpMatrix, x + 1, y + 1, p3, zoomFactor); tmpP4 = this .FindOrCreateNewPoint(tmpMatrix, x, y + 1, p4, zoomFactor); edge1 = this .FindOrCreateNewEdge(tmpEdges, tmpMatrix, x, y + 1, x + 1, y + 1); edge2 = this .FindOrCreateNewEdge(tmpEdges, tmpMatrix, x, y, x, y + 1); edge3 = this .FindOrCreateNewEdge(tmpEdges, tmpMatrix, x, y, x + 1, y); edge4 = this .FindOrCreateNewEdge(tmpEdges, tmpMatrix, x + 1, y, x + 1, y + 1); edge5 = this .FindOrCreateNewEdge(tmpEdges, tmpMatrix, x, y, x + 1, y + 1); triangle1.A = edge1; triangle1.B = edge2; triangle1.C = edge5; triangle2.A = edge3; triangle2.B = edge4; triangle2.C = edge5; result.Add(triangle1); result.Add(triangle2); } } matrix是原始数据,在我的项目是照度点矩阵,step的作用是由原始数据中生成三角网时可以跳过step个点取一点。 有了三角网,接下来应该去进行等值线追踪了 有了三角网,接下来应该去进行等值线追踪了。 在我一开始的想法中,绘制等值线肯定就是把所有具有相同值的点连接起来就OK了,一想那不是一个蜘蛛网嘛,不同高程值的等值线都交叉了,那还叫什么等值线?了解到使用三角网剖分方法来生成等值线后,又是以为三角形延着三个点来游走就能得到等值线,又一想不行,因为如果值恰好在某一个点上时,那到底向哪条边游走呢?最重要的是,如果我要绘制照度为500的等值线,如果三角网中所有的点上的值没有500怎么办?那岂不是很滑稽:有600,有400,但是没有500这条线,这就好比有爷爷,也有孙子,但是没有儿子,那孙子是从哪儿来的呢?这才发现原来人家早提出方案了,使用插值法,更重要一点是论文中还特意使用等值线不从三角网上的点通过,为什么?因为从点上过的话,等值线的方向不好计算。 知道了等值线都是与三角形的边交叉后(当然,除了三角网的边框),这里还需要对等值线的另一个特性有个了解,即等值线按照其类型可以分为开放式和封闭式。所谓开放式是指,等值线的两个端点都在三角网的边框上,所谓封闭式是指,等值线是一个闭环,起点同时也是终点。 假定我们需要绘制高程值为500的等值线,那么在三角网中追踪等值线的基本思路是:先查找所有开放式等值线,然后再查找封闭式等值线。开放式等值线的具体作法为: 1. 将所有至少两个边都包含500这个值的三角形找出来放到一个列表valTrangles中,判断一个边上是否包含500这个值的点的方法: 图1 ? //为了避免边的两个端点的值正好等于指定值,从而造成难以判断等值线方向 //因此当端点值正好等于指定值时给端点值减去一个很小的正数,使指定值所 //在的点和两个端点不重合 decimal tolerance = 0.000001M; decimal value1 = edge.P1.Value, value2 = edge.P2.Value; if (value1 == value) { value1 -= tolerance; } if (value2 == value) { value2 -= tolerance; } return (value1 - value) * (value2 - value) = 0; 从图1和上面的代码中我们可以清楚的看到,如果V点介于V1和V2的话,那么value1 - value和value2 - value中必然有一个是负数,从而使代码返回值是true。 2. 从valTrangle中任意一个三角网边框边的三角形开始,判断这个边上是否具有500这个点,如果有,记录当前三角形以及当前边,并将该点添加到等值线构成点列表(var Contour = new ListVPoint())中,开始在这个三角形剩下的两个边上查找“出路”,(我的代码没有处理剩余的两条边上同时具有指定高程值的情况,但是直至目前没有出现不合理的等值线。)。这里有一个问题:知道一个边上有某个点,那么这个点的坐标是多少呢?因为我们的三角形的边长度比较小,所以我们可以利用两个点,通过线性关系得到指定值所在的坐标: ? var factor = ( float )((value - edge.P1.Value) / (edge.P2.Value - edge.P1.Value)); result.X = edge.P1.X + factor * (edge.P2.X - edge.P1.X); result.Y = edge.P1.Y + factor * (edge.P2.Y - edge.P1.Y); 3. 找到第一个“出口”点后,将该点添加到Contour中,并将当前三角形标记为-1,即临时已使用。然后使用当前边为“出口”点所在的边,使当前三角形为当前边的另一个所属三角形,还记得我们Edge类型中的Trangle1和Trangle2的定义吧。 ? currentTrangle = currentTrangle == currentEdge.Triangle1 ? currentEdge.Triangle2 : currentEdge.Triangle1; 4.在得到当前三角形,当前边后循环调用第3步,又会得到第3个等值线的点,依照这个方法,最终会有两种可能,a:某个边只有属于一个三角形(到达三角网的边框了);b:到达的三角形的UsedStatus为-1(这个三角形已经被前面的步骤使用过了)。这时退出循环。 这个时候,对于某个值的开放式的等值线就已经查找完了,这里使用一个边所属的两个三角形这种方法,可以在很大程度上降低使用点阵获取三角形和边的麻烦,那种方法虽然谈不上非常困难,但是一大堆的索引值足以使你心烦意乱。 查找封闭式等值的方法和开放式基本上没有区别,所不同的只是退出条件为找到的出口和等值线第一点的坐标相同而已。 这里给出核心代码(同时包含开放式和封闭式): ? private IListVPoint GetContour(IListTriangle triangles, Edge currentEdge, decimal value) { IListVPoint result = new ListVPoint(); var firstEdge = currentEdge; var currentTrangle = currentEdge.Triangle1; var currentPoint = this .TryGetContourPoints(currentEdge, value); result.Add(currentPoint); while ( true ) { currentTrangle.UsedStatus = -1; var edges = currentTrangle.Get2OtherEdge(currentEdge); currentPoint = this .TryGetContourPoints(edges , value); int useEdge = 0; if (currentPoint == null ) { currentPoint = this .TryGetContourPoints(edges , value); useEdge = 1; } result.Add(currentPoint); currentEdge = edges ; currentTrangle = currentTrangle == currentEdge.Triangle1 ? currentEdge.Triangle2 : currentEdge.Triangle1; if (currentTrangle == null || currentTrangle.UsedStatus == -1 || currentTrangle.UsedStatus == 1) { int usedStatus = 0; //开放式等值线 if ((firstEdge.IsBorder currentEdge.IsBorder result.Count = 2) || //封闭式等值线 (firstEdge.IsBorder == false currentTrangle != null currentTrangle.UsedStatus == -1 result.Count = 6 currentEdge == firstEdge)) { usedStatus = 1; this .CompleteContour(triangles, usedStatus); } else { this .CompleteContour(triangles, usedStatus); result = null ; } break ; } } return result; } 当查找结束后,如果等值线合法,即开放式到达了边框,封闭式到达了起点,则“提交”一下候选的三角形,使用从临时使用状态回退到未使用状态1,不再参与下次查找,反之,则“回滚”一下候选三角形,使用从临时使用状态回退到未使用状态0。这里还有一个“保险”,开放式的等值线,一定是大于等于2个点,封闭式的等值线,一定是大于等于6个点。 感谢《 等值线标注的一种算法探讨 》一文的作者,我正在是使用这篇论文中的重要算法指导了我的工作。 首先标注那些小的封闭式的等值线。 这里我也没有想出来好的方法,就使用方法,找出封闭式等值线中点坐标X最小值,Y最小值,X最大值及最大值;如果XMax - XMin 指定值以及YMax - YMin 指定值,则在P((XMin + XMax) / 2,(YMin + YMax) / 2)处把等值线的值画出来,如下代码所示: ? if (points .X == points .X points .Y == points .Y) { float xMin = points.MinVPoint(p = p.X), xMax = points.MaxVPoint(p = p.X), yMin = points.MinVPoint(p = p.Y), yMax = points.MaxVPoint(p = p.Y); if (xMax - xMin 25 yMax - yMin 25) { g.DrawString(value.ToString(), font, brush, (xMax + xMin - sf.Width) / 2f, (yMax + yMin - sf.Height) / 2f); continue ; } } 效果如图中红圈所示: 图1 把过于小的封闭式等值线标注后,大点的封闭式和开放式等值线的标注方法就一样了。 根据《等值线标注的一种算法探讨》一文中的指导思想,需要把曲线转换为折线或是多边形,啥意思咧,看下图所示: 如果我告诉你,现在如果每条线段的长度大于高程值的字符串(value.ToString())所需的长度 - 空白长度,就在那里画一个标注应该就问题不大了吧,对于那些太小的线段,就不要去画了,小到不像话的封闭式等值线,我们在上一步中已经处理过了。 那么如何把一条曲线转换为折线或是多边形的呢?《等值线标注的一种算法探讨》一文中告诉我们(见论文2.2):“设一个dif参数,用来控制多边形近似等值线的误差,将曲线上的点从第0个开始,偶数点相连作为多边形的边,依次查看等值线上的3个点,如果中间的等值点(奇数点)到多边形的边的距离小于参数dif,舍弃中间的等值点,否则保存中间的点。如此不断循环,直到最后生成的多边形的边数不在(再 think8848注)减少为止。这样就得到了等值线近似的多边形。”dif越大,表明曲线越陡,dif越小,表明曲线越平缓。 算法代码示例: ? int n = (points .X != points .X points .Y != points .Y) ? points.Count : points.Count - 1, minEdge = n + 1, k = n; float tolerance = 8f; var indexes = new List int (); for ( int i = 0; i points.Count; i++) { indexes.Add(i); } while (k minEdge) { minEdge = k; var p = 0; while (p minEdge - (n % 2 == 0 ? 3 : 2)) { float straight = this .GetPointToStraight(points ], points ], points ]); if (Math.Abs(straight) tolerance) { indexes = -1; k -= 1; } p++; p++; } indexes = indexes.Where int (index = index != -1).ToList int (); } 另附求点到直线的距离的方法,下列代表示例如何计算p点到直线p1p2最短距离 ? private float GetPointToStraight(VPoint p, VPoint p1, VPoint p2) { if (p1.X == p2.X) { return ( float )Math.Abs(p1.Y - p2.Y); } if (p1.Y == p2.Y) { return ( float )Math.Abs(p1.X - p2.X); } double k = (p2.Y - p1.Y) / (p2.X - p1.X); double c = (p2.X * p1.Y - p1.X * p2.Y) / (p2.X - p1.X); return ( float )((k * p.X - p.Y + c) / (Math.Sqrt(k * k + 1))); } 最终indexes里面保存了在等值点列表中构成多边形的点的索引。 在本文的最后,我们再来谈谈如果标注字符旋转的问题,有一条线段作为基准,将画布旋转与线段倾斜角度相同度数应该不是一件难事,是的,.NET很容易就能做到,唯一的问题是如何线段的倾斜角度: ? var alpha = ( float )(Math.Atan((p2.Y - p1.Y) / (p2.X - p1.X)) * 180 / Math.PI); 就是这个角度了,三角函数已经忘了的兄弟可以到网上重新温习下高中数学,呵呵。 ? g.TranslateTransform(xOffset, yOffset); g.RotateTransform(alpha); g.DrawString(value.ToString(), font, brush, new PointF(0, 0)); g.RotateTransform(-alpha); g.TranslateTransform(-xOffset, -yOffset); 变换坐标系的原点到将要绘制标注的左上角位置,然后旋转画布,(注意,角度为正时为顺时针,角度为负时为逆时针方向,好像和我们数学课中的方向不一致。)绘制完成后,再把坐标系归位,循环,直至将所有的标注都绘制完成。 至此,使用C#绘制等值线的工作基本完成。
个人分类: 技术|7276 次阅读|1 个评论
余志伟老师的点滴追忆
热度 6 yuxuetaoxp 2011-12-15 20:58
第一次听说余老师,可以从我研究等值线的插值算法说起。那是我在研一的时候,由于项目的需要,需要加入根据雨量站统计的雨量值绘制研究区等值线的功能,那时候,我把国内所有关于等值线的文章几乎都大致看了一遍,有些收获也有些感触,那时候我看到矿大学报 1987 年发表的一篇关于等值线绘制的文章,署名作者是余志伟,单位也是中国矿业大学,当时自己也感叹认为此人物应该说是在地学界应用等值线比较先驱的人物了。所以这个名字也就在自己的脑海里有了一个根。 有次我跟我硕士导师奚砚涛奚老师聊天,说到了我正在研究的等值线。当时几乎我的脑海里几乎就是这个算法的实现。正巧,那天奚老师跟我说:“你可以看看余志伟余老师的文章,他是我们系研究等值线比较早的,而且首次提出了曲面样条插值算法,而且他认为那是他的最得意之作”。奚老师还说到,“余老师人非常好,他跟谭老师(谭海樵,我硕士导师的导师,即我师爷)关系非常好,那年我留在矿大任教就是余老师让我留下来的,而且芮老师(芮小平)就是他很早的硕士学生”。奚老师还说,余老师这几年身体不大好,他 01 年去北京后 04 年还来过矿大,后面就再没有来过了,如果下次去北京,一定要去看望看望他。奚老师这么一说,让我跟余老师的距离几乎从万米缩小到了几米的距离。 就这样,我也慢慢了解到了余老师的一点点的个人信息。到 10 年我硕士将要毕业准备考博的时候,机缘巧合,芮老师通过余老师将我推荐到了中国气象局程明虎老师那边,由于他们考数字信号处理、概率论和随机过程等课程,其中的数字信号处理和随机过程两门课程我都没有学过,很遗憾最终没有录取。由于我知道考得很糟,而其他院校的考博报名时间即将截止,最后我选择了考取余老师这边的博士。 那天第一次见余老师,是一天的上午,因为那天余老师一二节有课,故我在十点后才去余老师的办公室综合楼 520 室。跟余老师虽从未谋过面,但是一直感觉跟余老师没有任何的距离,就像非常要好的朋友虽久未谋面,但还是那样的肝胆相照一样。那天跟余老师聊得很嗨,说了一些我的个人情况,包括我写过什么论文和做过什么项目什么的。最后余老师请我吃饭时还聊到了当前 GIS 发展的瓶颈,我说当前 GIS 在煤矿上根本用不上,比如说算储量吧,顶多算算后就再没 GIS 的什么事了。余老师马上反驳我说,目前 GIS 的应用情况为什么是这样,那是因为我们没有把 GIS 做成动态或者时态的,如果随着煤矿的开采,我们能够实时地观测到地下地质条件和水文状况,那将是目前煤矿上所亟需的,而不是我们随便画画等值线和做个储量计算所能解决问题的。听余老师的这番话,让我感觉 GIS 的发展前景还是蛮大的,只是我们目前做得还不够完善。 就这样,点点滴滴,让我渐渐地清晰地了解余老师,而且我也顺利考上了余老师的博士。在我硕士答辩的时候,在答辩的前几天我听说了余老师要来,心中暗想,余老师这次来矿大,我觉得其中的一个重要的原因就是因为我要读他的博士。在我们的硕士论文答辩中,我还记得还有一个小插曲,我的一个同学做论文做得不是很好,被余老师批了。那个同学的论文题目是基于全球尺度的碳通量的计算,而实际中他只采用了加拿大的数据,故余老师的意思是让他把题目改小点,题目太大,而那个同学根本不识趣,拿他的论文指导老师陈报章老师来压所有的答辩评委老师,结果被余老师的一句你是让他答辩还是让我们答辩给顶了回去。从这件事中也印证了我的一个对余老师的看法,即对做学问追求一丝不苟。那天我的答辩比较顺利,老师问的问题我也基本都答上来了,毕竟我做这块做了三年,有了很好的基础,对整个行业也比较了解。我的论文也顺利地被评为了校级优秀硕士论文。那天我的答辩感言是“余老师身体一直都不太好,而他能够从北京过来参与我的硕士论文答辩,一定与我将要读他的博士有一定的关系,故我非常感激余老师能够过来参与我的硕士论文答辩。”虽然余老师只是笑着抽了口烟说:“不是这么回事,跟你没关系”。可我知道余老师一向就是一个比较倔脾气比较怪的人,口上虽不说,但是心里有这层的意思。 到后来我来了北京,余老师也经常找我谈话,说说目前的研究进度。我也经常找他谈论我的一些想法,一些研究的 idea 。偶尔也会向他请教问题。我一直说我要写一篇大作,发一篇超级牛的文章,余老师也默默地支持着我的想法。就在余老师离开我们的前几天,我还向他请教关于立方体插值的算法,他把他认为我会用到的文章全部毫无保留地拷给了我,让我多看。目前我已经想到了一个新的插值算法,可我想再向你请教时,你已离我们而去了。 最后一次见余老师,是在 2011 年 12 月 8 日的晚上,那天我注册参加的一个会议要有两个人去,我想问问余老师会安排谁陪我去,我说了杜晓敏师兄,余老师说“杜师兄应该没时间,我陪你去吧。我还没有去过哈尔滨呢, 1 月份的话哈尔滨还是比较冷的。行了,到时候我陪你去,你不用操心了。”就这几句,竟成了余老师留给我的最终遗言,您说的一起去哈尔滨的愿望也只能有我去带您实现了。 余老师虽离我们远去了,但是您探索科学的精神一直在激励着我,您也一定在天堂看着您的学生成长,成为您的衣钵传人。愿余老师在天堂安好,愿余老师家人能够从痛苦中走出来,能够去实现余老师的遗愿,让余老师安息、安眠。 您的博 11 级学生:于雪涛
个人分类: 生活点滴|4913 次阅读|12 个评论
C#实现基于三角网的等值线追踪及填充算法
热度 15 moustudio 2010-3-16 15:15
Captain Dialog 2010-03 1 引言介绍 等值线是一种离散数据的图形表示方法,在水利、土木、地质、石油勘探等工程和技术领域内广泛的应用。常规的等值线绘制通常采用网格法,其绘制的步骤一般为:离散数据网 格化;等值点的计算;等值线的追踪;光滑和标记等值线等。等值线图的显示方式一般有两种: (1) 等值线显示,即采用线条上加注数值标记的方式显示数据,这种方式的特点是简捷; (2) 采用彩色填充的方法来显示数据,既用不同的颜色来显示不同的数据,这种方法的特点是比较直观。两种方法的计算机实现也各不相同,一般来说,它们都需要将数据进行要用 到的网格网格化。第一种方法必须进行等值线的追踪、光滑和标记等值线。而第二种方法可以在追踪出等值线的基础上进行,也可以不做等值线的追踪直接在网格数据上进行操作。 方法实现的难易程度各不相同。 本文参考了文献【 1 】填充等值线的方法原理,并针对文中的寻找出来用于填充的多边形不是最小多边形进行改进。基于计算机图形学的基于左转算法的多边形自动搜索方法可以有效的识别最简单多边形,并结合方位角计算排除重复的多边形。借此提高了填充的效率。并且,结合计算机绘图的原理,本文还在填充的基础上,实现了等值线的追踪,相对经典的等值线追踪算法,此法更加简单,易于编程。 针对非规则数据(即散乱点数据),本文同样提出了构造等值线的方法,首先采用计算机图形学中的二维。。。生成非规则的三角形网格,然后针对每个三角形使用类似规则矩形的处理方式,可顺利解决不规则等值线的绘制问题。 2 散乱点数据的三角网生成 针对不规则形状的等值面,需要对特殊情况进行边界的处理,用规则网格情况下,经常进行等值线的边界进行先扩充再裁剪等方式,结果使得等值面的边界为锯齿状边界;然后利用三角形来构造边界,则会使得边界区域光滑,从而消除锯齿的产生。因此,在处理由散乱点数据构成的不规则等值面时,需要首先对所有的数据点进行三角网的生成,用三角网来构建等值面。 3设计的技术细节关键点 1 、等值线的标注: 标注的数值通常写在曲线比较平坦的部分,方法如下: (1) 寻找标注位置。为寻找曲线比较平坦的部分,依次找出 3 点,计算 3 点间两线段的夹角,若夹角大于 120~ ,则认为该处适合标注。 (2) 调整等值点顺序。要求标注的数值写至等值线中间时,需断开原等值线。若原等值线是非闭合曲线,则被分成两段,原来的等值线起始点和终止点不变,但在切断处增加一个新的终止点和一个新的起始点;若原等值线是闭合曲线,则可当作非闭合曲线处理,把起始点和终止点位置调整到切断的位置。 2 、等值线的追踪: 从图形角度分析,可利用当前点的颜色值,确定出一条等值线。 3 、等值线的光滑处理: 三次 B 样条函数处理插值处理 等值线的计算 判断格网的一条边是否与雨量值为 的等高线相交,要看这条边的两个端点的雨量值是否含有这个 值,例如点 A 、 B 是某个三角网的两个顶点,其雨量值分别为 及 ,则边 AB 是否与雨量值为 的等高线相交,应判断不等式 ( 一 )(He 一 ) 0 是否满足,就可以知道边 AB 是否与等高线相交。以往的此类算法为 J 编程的简单,对于 W= 或 : He 的情况, 采用将雨量值增加微小值的方法,对于此类情况作了专门处理,不需要改动雨量值而是直接追踪,使生成的雨量线更为准确。如图 1 ,穿过点 A 的雨量线可以看成如虚线的雨量线来处理,只是点 F 和点 G 退化到了点 A ,这样追踪时就要经过边 AB 、边 AE ,相当于绕着顶点转到了边 DE 。 判断一条边 AB 是否与等值线 W 相交的依据为: 升级处理:计算机中,除法容易导致崩溃,因此采用乘法处理: 4 、以颜色表为等值线值 5 、程序实现 应该用类似于计算机图形学的东西来实现,并且针对三角网等图元,采用面向对象的思路,使得每个图形元素的属性值可以得到灵活的操作性。 点 --- 线 --- 三角形 相同属性的点构成线; 相同属性的线封闭构成块; 点又可分为普通的数据点(控制点)、等值点(用于构成等值线)。 线于是也可分成边界控制线(由普通电构成)和等值线(由等值点构成)。 功能设计: 点: Cmou_Point ID 号? 坐标属性( x , y )、具体的属性值( V ) 所属的线号( list_ID_Line )? 点的类型(原始数据点,等值点) Note :在设计过程中,为了保证数据的完整性,基于点构造出来的段、线、三角形等包含的点信息,一律采用点的 ID 号表示,这样,针对点的修改就可以很容易实现,修改一处,用到该点的其他地方就不用重复修改了。同样,线号、段号都是这个用处! 段: Cmou_Segment ID 号 包含的起点和终点 StartPoint 、 EndPoint 所属的线号 (list_ID_Line) 所属的三角形编号( list_ID_Triangle ) 段的类型(想对于三角形来说:普通边、等值边,即含有等值点的边) 线: Cmou_Line ID 号( ID ) 线的类型(三角形边线,等值线,) 包含的数据点集( list_Point ) 线所属的三角形 ID 号 (list_ID_Triangle) 三角形: Cmou_Triangle ID 号 // 三角形包含的边线( list_Line ) 三角形包含的线段 : ( Cmou_Segment 类型) Edge0 、 Edge1 、 Edge2 利用线段构成了三角形 在三角形边上找到等值点后,三角形需要细分为多个三角形? 点集: list_Point_Scatter 原始的散乱点数据集 List_Point 所有用来生成三角形网格的数据点(包括原始数据点和新增的等值点) 三角形的等值线追踪时,可以通过三角形边线上的 ID_Triangle 追踪到与之相邻的三角形(两个相邻的三角形共边) 等值线?? Cmou_ContourLine ID 编号 关于三角形的追踪,可以构造出一条条的 Cmou_Line ,然后绘制显示即可。 // 关于三角形的填充,可以最后根据三角形所含的数据点的 V 值,构造出人意形状的多边形, // 然后按照等值线的数值进行填充。 关于等值线的填充,完全可以按照等值线的值的大小排列,由两条等值线构成多边形,然后再用对应的颜色填充。 多边形: Cmou_Polygon ID List_Point 包含的点集 Vstart 、 Vend 属性值信息 等值线追踪算法描述: (1) 非闭合等值线对应等值边系列的搜索。 ① 以某等值线的数值为依据,先在 边界边 上找到一条等值边作为起始边,若找不到这样的起始边,则说明该等值线闭合,执行 (2) 。将该边记录到 等值边系列 后,从 待搜索边系列 中剔除。 ② 在起始边相关三角形的相关边上寻找第二条等值边,找到后加入等值边系列,并将其从待搜索边系列中剔除。 ③从 待搜索三角形系列 中剔除等值线已经过的三角形,则最新找到的等值边只有一个相关三角形,在该三角形内搜索下一条等值边,找到后加入等值边系列,并将其从待搜索边系列中剔除。 ④递归执行③ , 直到找到的最后一条等值边是边界边 。 ⑤ 一条非闭合等值线的等值边系列搜索完毕。 (2) 闭合等值线对应等值边系列的搜索。 ① 任找一条 等值边 作为起始边,若找不到这样的起始边,则说明分析不出该等值线,应退出该等值线的分析而开始另一条等值线的搜索。将该边记录到等值边系列,并在该边上作一标志。 ②在起始边的一个相关三角形中寻找第二条等值边,找到后加入等值边系列,并将其从待搜索边系列中剔除。 ③从待搜索三角形系列中剔除等值线已经过的三角形,则最新找到的等值边只有一个相关三角形,在该三角形内搜索下一条等值边,找到后加入等值边系列,并将其从待搜索边系列中剔除。 ④递归执行③ ,直到找到最后一条等值边是起始边。 ⑤一条闭合等值线的等值边系列搜索完毕。 (3) 循环执行 (1) 、 (2) 直到找到所有等值线对应的等值边系列。用线性插值的方法求等值边上对应等值点的坐标。将这些等值点有序地连起来即实现了等值线的自动分析。 等值线追踪相关【 3 】: 对每条等值线值 ,等值线的游走可按如下 3 个步骤进行。 (1) 搜索存在等值点的所有三角形:首先在三角剖分后得到的三角形链表中进行搜索,找出所有存在等值点的三角形,并将此三角形添加到一个临时建立的三角形链表中供以后使用 ( 此链表记为 TriCheckedAngleList) 。 (2) 搜索开曲线等值线:从凸边三角形链表 (TriChimbAngleList) 中搜索等值点存在的三角形,找出等值线的起始点,记下点的位置和点所在的边 ( 即两个三角形的公共边 ) 。然后以此 公共边 在 TfiCheckedAngleList 链表中寻找包含此边的新 三角形 ,找到后,在新三角形中进行等值点的判断和计算,记录等值点和新的公共边,并以新公共边按上述方法继续在 TriCheckedAngleList 链表中游走,直至找不到包含新公共边的三角形为止。需要注意的 是在游走过程中,必须将每次找到的新三角形从 TfiChecke3dAngleList 中删除,否则等值点会出现重复。按此方法进行等值线游走后,就得到了两端经过三角形网格边界的开曲线等值线。 (3) 搜索闭合等值线:完成开曲线等值线搜索后,在剩下的 TriCheckedAngleList 链表的三角形中取第一个三角形,在其边上找出等值点后,按步骤 (2) 同样的方法进行等值线游走,直到 TriCheekedAngleList 链表的三角形为空,这样就可搜索出 值的所有闭合等值线。 Note :关键在于如何处理三角形中搜索等值边的递归关系。 参考文献 : 庞世明 , 蔡玉华 , 靳文芳 . 等值线图的彩色填充方法 . 计算机应用 . 2004.1 Vol.24,NO.1, 60~62. 李顺新 , 刘 俊 , 陈建勋等 . 利用不规则三角网法绘制三峡雨量等值线图 . 2005.4,Vol.36,No.4, 44~46. 赵 伟 , 赵卓宁 , 李五生 . 一种有效的离散数据场等值线生成方法 . 成都信息工程学院学报 .2007.2,Vol.22,NO.1 115~121. 2010-2-10 等值线的追踪工作基本完成 2010-2-11 关于等值线的填充: 需要制定基于单个三角形的填充,将所有的数据点(包括原始的散乱点、等值点、以及将要添加进去的圆滑插值点)进行任意形状的填充。 由此,形成有三套数据点集。每个三角形进行填充的时候都是需要查询这些点集,找出自己的点,然后根据点的属性值构造出任意多边形,再进行填充。 终于找到错误的原因了!!! 初步确定问题为三角形中的等值点没有连续。如:上面的 44 号三角形( 102-106-107 )中顶点和等值点是: 93.877 , 94.812 , 94.812 , 96.260 , 98.675 (除了红色为等值点,剩下的为顶点) 而控制它们的等值线值为: 91.682 , 94.812 , 97.943 于是,最后的一个点 98.675 成了孤立点,不能参与到三角形的填充搜索中。 解决办法:仔细调试后,发现问题应该是出在等值线的追踪过程中。因此,修改了追踪算法中的一个条件,使得在一般情况下,不会出现上述问题。 但是,上面的那种情况没有做处理,只是认为正常的等值线追踪情况下,不可能出现那样的等值点不连续的情况。 关于等值线值相应的颜色填充 相应的颜色填充,可根据前面等值线填充中的判断等值点的值属于受哪个等值线值控制,取对应的颜色。下标范围是 [ 最小,最大 ) 下一步工作: 1、 等值线的光滑 2、 等值线的标注 3、 等值线的颜色改进 4、 任意形状处理 2010-2-23 改进之增加颜色选择按钮 修改等值线颜色 修改等值线值的字体颜色 修改等值线的背景颜色(在等值线填充的情况下无效) 增加了图片保存功能
个人分类: 算法研究|26133 次阅读|9 个评论

Archiver|手机版|科学网 ( 京ICP备07017567号-12 )

GMT+8, 2024-5-12 08:27

Powered by ScienceNet.cn

Copyright © 2007- 中国科学报社

返回顶部