十年匠心定制 · 商业建站与技术教学双线并行 咨询热线:400-886-1026 service@lmnt.cn
ARTICLE DETAIL

资讯详情

深耕网站建设与运营推广的一线实战洞察。

Excel与KML互转及CZML动画生成:无人机航迹可视化实践指南

Excel与KML互转及CZML动画生成:无人机航迹可视化实践指南 简介这是一份GIS专业工具软件包面向地理测绘、无人机巡检、规划展示等场景解决KML与Excel数据批量互转、无人机照片与传感器数据处理、以及KML/CZML时间轴动画生成等常见难题。压缩包共742个文件约66.75MB核心包含可执行的exe主程序、支撑运行的大量dll动态库、用于坐标与投影参数管理的csv数据文件以及kml、dxf、dae等地理模型样例同时附带数百张png界面截图和txt说明文档方便用户对照查看操作流程。目前已有736人学习下载。借助该包可快速部署“大奇GIS专家”工具将Excel表格一键转为地图标注也能将无人机采集的影像整理为数字高程模型或正射影像并导出带时间维度的CZML动画供Web端展示。随包附带的坐标系、投影参数和字段配置文件对理解GIS底层数据组织方式也很有参考价值。1. 大奇GIS专家到底在解决哪个工作流里的问题一个航次飞完机载POS导出一张Excel表两三百行经纬度、高度、姿态角老板要你把飞行轨迹放在Google Earth或Cesium里回放。这时候最尴尬的不是“没有工具”而是手头同时缺三样东西一张能直接吃表格坐标的KML一条能从无人机Excel数据里梳理出的航迹以及一套能把轨迹按时间轴播出来的CZML动画。大奇GIS专家这类工具包瞄准的正是这个链条——让Excel里的坐标“活”到地图上再让地图上的轨迹“动”起来。面向的人群很明确GIS开发工程师、无人机数据处理岗位、三维可视化实施的人。凡是手边躺着“坐标系一堆表格”的人都会遇到这个标题里描述的三个中转站。下面按我平时处理数据的顺序把每一段拉通。2. KML与Excel相互转换先把表格变成地图再把地图变回表格我平时处理Excel与KML互转时不太会把“大奇GIS专家”当成黑盒而是直接按下面这段逻辑自己生成中间文件。KML本质上是一个以Placemark为单位的XML文档一个Excel行对应一个Placemark列与节点一一映射。“相互转换”这件事用不着重型GIS桌面软件一段Python脚本就能闭环。关键要处理三件事坐标字段的抽取、几何类型的映射、以及KML里那两个国际通行问题坐标系顺序和UTF-8编码。2.1 从Excel到KML的最小生成模板gis坐标转点先用最简单的“点表”说起这是GIS教程里最常出现的例子。Excel里至少要有lon、lat两列最好再有一列作为Placemark的name。用openpyxl读文件用字符串拼接把每个点写成Placemark。import openpyxl from xml.sax.saxutils import escape def excel_points_to_kml(xlsx_path, sheet, lon_col, lat_col, name_col, out_kml): wb openpyxl.load_workbook(xlsx_path, data_onlyTrue) ws wb[sheet] kml [] kml.append(?xml version1.0 encodingUTF-8?) kml.append(kml xmlnshttp://www.opengis.net/kml/2.2) kml.append(Documentnameexcel_points/name) for row in ws.iter_rows(min_row2, values_onlyTrue): lon row[lon_col] lat row[lat_col] if lon is None or lat is None: continue name escape(str(row[name_col])) if name_col is not None else # KML坐标顺序是 经度,纬度,高度 kml.append(fPlacemarkname{name}/name fPointcoordinates{lon:.8f},{lat:.8f},0/coordinates/Point f/Placemark) kml.append(/Document/kml) with open(out_kml, w, encodingutf-8) as fp: fp.write(\n.join(kml)) # 参数第0列lon第1列lat第2列为name excel_points_to_kml(pos.xlsx, Sheet1, 0, 1, 2, points.kml)这段代码的参数选择我一般这样定lon_col、lat_col直接用列的索引而不是列名避免不同版本Excel表头不统一用data_onlyTrue是为了防止公式列读出来是None。输出用lon:.8f度单位下8位小数约等于1厘米足够无人机航点使用如果只是宏观轨迹6位小数就够。生产环境里数据量大时我不会用字符串拼接而是直接拿geopandas.GeoDataFrame.to_file(driverKML)做批处理但小体积的Excel转点上面这段没有第三方GIS依赖复制过去就能跑。如果要做批量出图把这个函数放进for sheet in wb.sheetnames循环里每张表导出一个KML就是最简单的批量出图方式。2.2 反方向解析把KML的Placemark坐标抽回Excel反向转换要处理的除了Point还有LineString和Polygon。用xml.etree.ElementTree遍历命名空间里的Placemark把坐标字符串按“空格分点、逗号分经纬度”切回来写到openpyxl的单元格里。import openpyxl, xml.etree.ElementTree as ET NS {kml: http://www.opengis.net/kml/2.2} def kml_coords_to_excel(kml_path, xlsx_path, sheetSheet1): root ET.parse(kml_path).getroot() wb openpyxl.Workbook() ws wb.active ws.title sheet ws.append([name, geom_type, lon, lat, height]) for pm in root.findall(.//kml:Placemark, NS): name_el pm.find(kml:name, NS) name name_el.text if name_el is not None else point pm.find(.//kml:Point/kml:coordinates, NS) line pm.find(.//kml:LineString/kml:coordinates, NS) poly pm.find(.//kml:Polygon/kml:outerBoundaryIs/kml:LinearRing/kml:coordinates, NS) coord_text None geom_type if point is not None: coord_text point.text geom_type Point elif line is not None: coord_text line.text geom_type LineString elif poly is not None: coord_text poly.text geom_type Polygon if coord_text: for part in coord_text.strip().split(): lon, lat, h part.split(,) ws.append([name, geom_type, float(lon), float(lat), float(h)]) kml_coords_to_excel(drone.kml, output.xlsx)逻辑说明ElementTree的findall里带冒号的标签必须通过命名空间字典NS来解析漏了这一步查找结果永远是空列表Polygon我只取了outerBoundaryIs如果遇到带洞的面需要再加innerBoundaryIs的遍历。输出时把name重复写到每一行是为了以后在Excel里用多条件筛选时还能按元素还原几何。如果要在Excel里按名称反查某个Placemark做法和“python查找excel中字符串”无异遍历ws.rows逐行比较即可。2.3 三个必调参数与三个高频坑先把参数表放这里后面所有脚本都会引用它。参数项推荐值说明坐标顺序longitude,latitude,heightKML规范是经度在前很多Excel表刚好相反高度单位相对WGS84椭球高单位米Google Earth默认显示海平面高度编码UTF-8换成GBK导出地图客户端打开会乱码小数位精度6~8位6位约0.1m8位约1cmXML声明version1.0 encodingUTF-8少了它部分浏览器直接报解析错误坑一坐标反了。KML里是“经度,纬度”而CAD图纸、全站仪导出经常是“纬度,经度”或“Y,X”。热搜里那个“cad到gis 6位坐标转换”是另一个经典坑CAD内部用的可能是毫米级投影坐标直接除1000变成米再用PROJ转地理坐标和KML的经纬度不是同一个坐标系别混用。坑二Excel无法粘贴数据。如果你试图把几千行坐标从一张表复制到另一张大概率看到Excel转圈后没有任何反应。这通常不是文件坏了而是剪贴板里格式太多。热搜里大量“excel无法复制粘贴”“excel复制粘贴不了”的解法都一样用openpyxl直接写单元格而不是经过系统剪贴板。上面两段脚本走的就是这条不经过剪贴板的路。坑三时间字段在KML里不会被Excel自动识别。when、begin、end这些时间字符串写回Excel后默认变成文本做时间筛选时要先分列或DATEVALUE转换。这个问题到第四章的动画生成里会更加明显。3. 无人机数据处理把POS记录变成可编辑的KML航迹无人机数据处理在GIS开发里通常不是“正射影像拼接”那种重计算环节而是先把POS表变干净再变成地图软件能用的航迹。POSPosition and Orientation System记录的核心字段就七个照片名、经度、纬度、高度、航向角yaw、俯仰角pitch、翻滚角roll。这里只关心位置姿态角留给后续场景重建。3.1 一张典型的无人机POS表长什么样下表是DJI类飞控导出的常见结构不同厂家列名不同但信息同构。Excel字段含义常见单位处理注意filename影像文件名文本会出现在KML的name里lontitude经度度这个列名经常拼错但DJI导出就是这么写的latitude纬度度有的飞控先纬度后经度altitude相对起飞点高度米和绝对海拔高程有差值yaw航向角度转CZML的orientation时需要转弧度pitch/roll俯仰/翻滚度转CZML时可以暂时忽略datetime拍照时刻yyyy-MM-dd HH:mm:ss本地时间转CZML时统一转UTC把这张表和上一章的KML对应起来无人机数据处理的第一步其实是“表结构清洗”统一列名、统一经度列在前、把时间字符串补齐格式。这也是为什么大奇GIS专家这类工具会把Excel转换和无人机数据放在同一个包里——两者共用的是同一套表格清洗逻辑。3.2 用Python把POS转成KML LineString这里给出读取CSV直接生成航迹的脚本。飞控导出的CSV通常不是标准UTF-8我加了一个encoding兜底坐标转换方面如果你的POS是RTK输出的WGS84直接用如果是大疆消费级导出部分第三方App输出的是GCJ02火星坐标需要先套一个转WGS84的偏移函数。import csv import math def gcj02_to_wgs84(lng, lat): # 一种常见的近似反向转换精度约1~3米够航迹展示用 a 6378245.0 ee 0.006693421622965943 dlat transform_lat(lng - 105.0, lat - 35.0) dlng transform_lng(lng - 105.0, lat - 35.0) radlat lat / 180.0 * math.pi magic math.sin(radlat) magic 1 - ee * magic * magic sqrtmagic math.sqrt(magic) dlat (dlat * 180.0) / ((a * (1 - ee)) / (magic * sqrtmagic) * math.pi) dlng (dlng * 180.0) / (a / sqrtmagic * math.cos(radlat) * math.pi) return lng - dlng, lat - dlat def transform_lat(x, y): # 实际项目里用官方接口或本地偏移表这里只做示意 return -100.0 2.0 * x 3.0 * y 0.2 * y * y 0.1 * x * y 0.2 * math.sqrt(abs(x)) def transform_lng(x, y): return 300.0 x 2.0 * y 0.1 * x * x 0.1 * x * y 0.1 * math.sqrt(abs(x)) def pos_to_kml(csv_path, out_kml): with open(csv_path, r, encodingutf-8-sig) as fp: rows list(csv.DictReader(fp)) coords [] for r in rows: lon, lat float(r[lon]), float(r[lat]) lon, lat gcj02_to_wgs84(lon, lat) # 确认是WGS84就删掉这行 height float(r.get(alt, 0) or 0) coords.append(f{lon:.8f},{lat:.8f},{height:.1f}) kml [?xml version1.0 encodingUTF-8?, kml xmlnshttp://www.opengis.net/kml/2.2Document, Placemarknamedrone_track/nameLineStringtessellate1/tessellate, fcoordinates{ .join(coords)}/coordinates, /LineString/Placemark/Document/kml] with open(out_kml, w, encodingutf-8) as f: f.write(\n.join(kml))逻辑说明tessellate1告诉Google Earth把线段贴合到地表否则直线会直接穿过地形utf-8-sig用于处理BOM头很多Windows导出的CSV不带BOM反而让中文字段名出问题高度这里直接取相对高度因为Google Earth显示的是海拔飞控的相对高度和绝对高程能差几十米。如果要拿这条轨迹去CZML里做高度曲线必须先做高程基准换算。3.3 航迹抽稀与平滑先设好两个阈值POS原始数据几百到几千行直接生成KML没问题但丢进Cesium做CZML动画时点数太多会拖慢插值。常见做法是先抽稀再平滑。抽稀我用最基础的“距离阈值”法相邻两点之间球面距离小于2米就跳过大于2米才保留这样既保留转弯特征又减少冗余点。def haversine(lon1, lat1, lon2, lat2): R 6371000 p1, p2 math.radians(lat1), math.radians(lat2) dp math.radians(lat2 - lat1) dl math.radians(lon2 - lon1) a math.sin(dp/2)**2 math.cos(p1)*math.cos(p2)*math.sin(dl/2)**2 return 2 * R * math.asin(math.sqrt(a)) def simplify(points, min_distance2.0): kept [points[0]] for p in points[1:]: if haversine(kept[-1][0], kept[-1][1], p[0], p[1]) min_distance: kept.append(p) return kept这个距离阈值是参数设计的核心2米在城区航线会留下大量点在空旷农田可以放到5米。我一般先跑一次统计CSV里的点间距分布再决定阈值。平滑则用三点滑动平均但注意别对经纬度分别平均然后放回——那样会在跨180度经线时产生瞬时偏移。3.4 从航点到航区GIS坐标成面与Excel多条件筛选除了航迹线无人机数据处理里还经常要出“作业范围”——把边界航点连成一个面。“gis坐标成面”就是这个需求。直接用shapely的convex_hull把POS经纬度点包成凸多边形再输出成KML Polygonfrom shapely.geometry import Point from shapely.ops import unary_union def points_to_boundary_kml(points): hull unary_union([Point(lon, lat) for lon, lat in points]).convex_hull ring list(hull.exterior.coords) coord_str .join(f{lon:.8f},{lat:.8f},0 for lon, lat in ring) return coord_str在这个阶段Excel多条件筛选也常被用来快速剔除异常点。比如用“航高 50 且 yaw 在 0 到 90”这种组合筛选把拐弯段或异常姿态段挑出来复查。这就是为什么我在2.2里特意把name重复写入每一行——回到Excel里可以对同一个Placemark下的所有坐标做条件格式和筛选不用再写SQL。4. KML和CZML动画生成从静态轨迹到Cesium飞行视角把KML轨迹变成CZML动画核心不是格式转换而是理解两者各自的时间模型。KML的gx:Track把时间和坐标分别放到when和gx:coord列表里CZML则用一个epoch起点加一组“相对秒数值”的数组来表达同一件事。两者的映射关系非常直接所以做转换时最值得花时间的是处理时区和插值参数。4.1 KML的gx:Track动画模型一段最小可用的gx:Track在Google Earth里可以播放点运动Placemark gx:Track when2024-05-01T08:00:00Z/when when2024-05-01T08:00:01Z/when gx:coord120.15 30.28 100/gx:coord gx:coord120.16 30.29 120/gx:coord /gx:Track /Placemark这里when必须按时间先后排序而且要带Z表示UTCgx:coord是“经度 纬度 高度”空格分隔不是逗号。很多工具生成的gx:Track时间没排序Google Earth会直接忽略整个轨道。写成CZML时这段数据的两个when会成为动画窗口的start和stop。4.2 用一段Python生成带时间的CZMLCZML的每个实体是JSON对象位置字段可以用cartographicDegrees直接写经纬度和高度配合epoch形成关键帧。下面这段代码把上一章抽稀后的航点转成CZMLimport json def pos_to_czml(points, start_time_utc2024-05-01T08:00:00Z): # points: [(lon, lat, height, elapsed_seconds)] epoch start_time_utc carto [] for lon, lat, h, t in points: carto.extend([t, lon, lat, h]) # 时间在前Cesium要求这个顺序 packet { id: drone, position: { epoch: epoch, cartographicDegrees: carto }, point: { color: {rgba: [255, 140, 0, 255]}, pixelSize: 8 }, orientation: { velocityReference: #drone.position } } with open(drone.czml, w, encodingutf-8) as f: json.dump([packet], f) # 假设每1秒一个点 pts [(120.15, 30.28, 100, i*1.0) for i in range(10)] pos_to_czml(pts)epoch之后数组的排列规则是“时间,经度,纬度,高度”这与KML的gx:coord顺序不同。velocityReference让模型自动朝向运动方向是做飞行视角动画的关键参数如果不需要朝向删掉orientation键就可以。interpolationAlgorithm默认是线性插值对无人机轨迹够用要做平滑曲线可以改成hermite但那样必须再给出tangent值否则Cesium会用默认值插出抖动的曲线。4.3 CZML加载时控制播放节奏的三个参数参数取值示例作用epoch2024-05-01T08:00:00Z锚定时间起点影响position里相对秒数availability2024-05-01T08:00:00Z/2024-05-01T08:02:00Z限定实体可见区间超出不显示interpolationAlgorithmlinear/hermite/lagrange位置、朝向在两个关键帧之间的插值方式shouldAnimateviewer.clock.shouldAnimate trueCesiumViewer里是否自动走时间轴在Cesium里加载CZML时CesiumViewer会自动扫描所有packet的availability并设置时间轴范围。如果发现动画没有按预期开始先检查epoch和start时间是否带时区后缀再检查间隔是否与数据点数量匹配。加载完成后记得给viewer.clock.shouldAnimate赋值true否则画面会停在第0秒。4.4 KML的when时间戳转CZML的偏移陷阱KML里时间通常是UTC的ISO8601字符串CZML同样用UTC。真正会坑人的是国内飞控导出的本地时间北京时间比UTC快8小时如果直接把2024-05-01 08:00:00当作UTC传给CZML动画里的点在真实回放时会晚8小时才运动。处理方式是先把本地时间转成带时区的UTC字符串from datetime import datetime, timedelta local datetime.strptime(2024-05-01 08:00:00, %Y-%m-%d %H:%M:%S) utc local - timedelta(hours8) print(utc.strftime(%Y-%m-%dT%H:%M:%SZ))同样的问题也出现在Excel里Excel的日期序列号不携带时区信息用openpyxl读出来的datetime对象默认是裸时间。我会在代码里强制做一次astimezone(timezone.utc)避免写进CZML后和Cesium的viewer时间轴差8小时导致动画跳到错误的日历日。5. KML、Excel与CZML处理结果的三步验证法最后这一章给出的是我每次处理完KML、Excel和CZML之后必跑的验证命令以及一张三源坐标自查表。多花一分钟做静态校验能省掉在Google Earth和Cesium之间来回切换的时间。5.1 ogr2ogr快速验证KML是否可被标准GIS库读取GDAL/OGR是绕不开的校验工具。拿到生成的KML第一件事就是用ogr2ogr把它读一遍看能否无损地转换回GeoJSONogr2ogr -f GeoJSON check.json drone.kml如果这条命令报Error 1: ... Invalid coordinate大概率是坐标字段里有空字符串、NaN或者“经度,纬度”顺序被写反。再把GeoJSON打开数一下Feature数量和原始Excel行数对比能发现ElementTree漏掉Placemark的情况。CZML不用ogr2ogr我用Python的json模块校验结构顺手检查epoch之后的数组长度是否是4的整数倍因为位置项必须是“时间lonlatheight”四元组。5.2 CZML文件结构与时间轴自检脚本下面这段脚本把Excel、KML、CZML三种来源的坐标全部读进来统一输出成经度范围、纬度范围和坐标顺序标记。凡是输出里lon范围落在(-180,180)正常区间、lat范围同样落在(-180,180)的都要怀疑经纬度反了lon范围正常但lat超过90度的一定是列读错。def check_coords(data, label): lons [p[0] for p in data] lats [p[1] for p in data] print(f{label}: lon[{min(lons):.6f},{max(lons):.6f}] flat[{min(lats):.6f},{max(lats):.6f}]) if max(lats) 90 or min(lons) -180 or max(lons) 180: print( !! 坐标顺序疑似颠倒或单位错误)5.3 三源坐标一致性核对表检查项通过标准失败时先看哪坐标值域lon∈[-180,180]lat∈[-90,90]是否是“纬度,经度”顺序时间排序KML的when单调递增gx:Track要求全局时间有序时区一致性所有时间戳都带Z本地时间转UTC高程参照相对高度和绝对高程不要混用飞控参数里的alt基准属性完整性Placemark数量等于源表行数ElementTree命名空间解析失败这套方式我没有做成大型校验库因为转换脚本每次面对的数据结构都不一样一张参数表比一个万能脚本更可靠。把上面的命令和范围自检加进处理流程之后再进Cesium或Google Earth出现的问题基本就只剩视觉调优了。本文还有配套的精品资源点击获取
返回列表