CN104537606B - 地理坐标投影变化方法 - Google Patents

地理坐标投影变化方法 Download PDF

Info

Publication number
CN104537606B
CN104537606B CN201410750917.3A CN201410750917A CN104537606B CN 104537606 B CN104537606 B CN 104537606B CN 201410750917 A CN201410750917 A CN 201410750917A CN 104537606 B CN104537606 B CN 104537606B
Authority
CN
China
Prior art keywords
projective transformation
coordinate
points
transformation region
region
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Expired - Fee Related
Application number
CN201410750917.3A
Other languages
English (en)
Other versions
CN104537606A (zh
Inventor
严泰来
叶思菁
朱德海
李林
张超
姚晓闯
杜振博
刘帝佑
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
China Agricultural University
Original Assignee
China Agricultural University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by China Agricultural University filed Critical China Agricultural University
Priority to CN201410750917.3A priority Critical patent/CN104537606B/zh
Publication of CN104537606A publication Critical patent/CN104537606A/zh
Application granted granted Critical
Publication of CN104537606B publication Critical patent/CN104537606B/zh
Expired - Fee Related legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T3/00Geometric image transformations in the plane of the image
    • G06T3/06Topological mapping of higher dimensional structures onto lower dimensional surfaces
    • G06T3/073Transforming surfaces of revolution to planar images, e.g. cylindrical surfaces to planar images
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/30Subject of image; Context of image processing
    • G06T2207/30181Earth observation

Landscapes

  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • Theoretical Computer Science (AREA)
  • Image Processing (AREA)

Abstract

本发明公开了地理坐标投影变化方法,所述方法通过提取变换区域的顶点坐标;计算顶点在目标投影坐标系下的同名点坐标;建立投影变换线性方程组;根据所述投影变换线性方程组计算所述投影变换区域内要素的目标投影坐标系下的坐标。本发明的方法在很小的局部区域内将经线与纬线分别用直线表示,将球面计算转化为平面计算,可以提高地理坐标投影变换的计算速度,并且可以通过控制局部区域的大小实现误差控制,能够以较少的成本满足涉及大量数据变换并且对变换精度与变换速度同时有较高要求的实际应用需求。

Description

地理坐标投影变化方法
技术领域
本发明涉及空间信息技术领域,更具体涉及一种地理坐标投影变化方法。
背景技术
经纬度与平面坐标是地理信息系统GIS的两种基本数据坐标体系。由GPS采集得来的点位数据以及传感器获取的卫星影像数据往往是经纬度数据(WGS84),而在行业应用领域,包括土地资源管理、交通运输管理等领域所处理和使用的图件资料,无论是矢量格式还是网格格式,往往是基于平面坐标的XY数据。因此在GPS数据,遥感数据以及图件资料等多源空间数据实际的集成应用过程中,经纬度数据投影变换的速度与精度是影响应用效率与准确性的关键因素之一。
经纬度数据与平面坐标XY数据之间存在着十分复杂的函数关系。对于坐标基准不变情况下的投影变换,根据投影参数,往往具备相应的投影正反算公式,转换公式非常复杂,以高斯投影正算为例,其计算公式为
式中,B为待转换点的大地纬度;l为待转换点与中央子午线间的经差;t=tanB;X为由赤道至纬度B的子午线弧长;η2=e'2cos2B,其中e'为参考椭球第二偏心率;N为卯酉圈曲率半径。这类转换在涉及坐标基准变化时会变得更为繁琐,现阶段国内外广泛采用的转换方法是相似变换法,其过程是将经纬度数据转换到空间直角坐标系上,通过一组或多组不同参考基准下两套坐标值的公共点计算两套坐标系统之间的平移、缩放、旋转关系,进而实现参考基准的变换,再在此基础上通过一定的变换公式实现投影变换。
近些年来,国内外学者对于相似变换模型和求取坐标转换参数的方法进行了大量的研究工作。早期阶段,莫洛琴斯基提出一种三参数转换模型,但随着观测精度的大幅提高,三参数转换法的精度已经完全不能满足测绘要求。目前通常采用的相似变换模型有七个转换参数,包括三个平移参数,三个旋转参数和一个尺度参数,例如布尔莎模型、武测模型和莫洛斯基模型等。为了进一步提高精度,有学者提出了包含两组旋转参数的转换模型,还有学者提出了参数可变的坐标转换模型,旨在降低地球重力等物理因素引起的残差。此外,相关的研究还包括优化公共点选择与参数求解算法等等。
现有技术的地理坐标投影变换方法虽然能在转换精度上满足实际应用的需求,并且应用起来方便灵活,但是在执行大批量矢量数据转换应用中,由于存在球面到平面的变换过程,计算较为复杂,在计算效率上即便借助云计算技术,也难以满足实际要求,这便在一定程度上影响了某些涉及大量数据变换并且对变换精度与变换速度同时有较高要求的业务实施。
发明内容
(一)要解决的技术问题
本发明要解决的技术问题是如何降低地理坐标投影变化计算方法的复杂程度,提高大数据量转换的效率与精度。
(二)技术方案
为了解决上述技术问题,本发明提供了一种地理坐标投影变化方法,包括以下步骤:
S101、将地理坐标系下的投影区域规划为一个二维的投影变换区域,所述投影变换区域为能够覆盖所述投影区域的面积最小的几何图形;
S102、提取所述投影变换区域的至少三个顶点的坐标;
S103、计算步骤S102的顶点坐标在目标投影坐标系下的同名点坐标;
S104、建立投影变换线性方程组;
S105、根据所述顶点坐标以及相应的同名点坐标,计算所述投影变换线性方程组的变换系数的值;
S106、利用所述投影变换线性方程组计算所述投影变换区域内要素在目标投影坐标系下的坐标。
优选地,所述投影变换区域为矩形,提取其四个顶点的顶点坐标:A(Ba,La)、B(Bb,Lb)、C(Bc,Lc)、D(Bd,Ld),其中A点和B点同一条纬线上,C点和D点在同一条纬线上,A点和D点在同一条经线上,B点和C点在同一条经线上;计算所述顶点A、B、D在目标投影坐标系下的同名点坐标A’(xa,ya),B’(xb,yb),D’(xd,yd)。
优选地,所述投影变换线性方程组为:
y=yaB(B-Ba)+θL(L-La)
其中(B,L),(x,y)为待投影要素分别在地理坐标系与目标投影坐标系下的坐标;θL,θB分别为所述变换系数。
优选地,所述变换系数的求解方法具体为:
将点A、A’、B、B’、D、D’的坐标值带入所述投影变换线性方程组,建立所述变换系数的计算公式:
θL=(yb-ya)(Bd-Ba)-(yd-ya)(Bb-Ba)/(Bd-Ba)(Lb-La)-(Bb-Ba)(Ld-La)
θB=(yd-ya)(Lb-La)-(yb-ya)(Ld-La)/(Bd-Ba)(Lb-La)-(Bb-Ba)(Ld-La)
将等值关系Ba=Bb,La=Ld代入所述变换系数的计算公式,可得到简化的公式:
θL=(yb-ya)/(Lb-La)
θB=(yd-ya)/(Bd-Ba)
根据点A、A’、B、B’、D、D’坐标值计算所述变换系数。
一种地理坐标投影变化方法,包括以下步骤:
S201、将投影区域划分为两个或两个以上的二维投影变换区域;
S202、为每个所述投影变换区域设置编码;
S203、对于每一个所述投影变换区域,提取其至少三个顶点的坐标;
S204、对于每一个所述投影变换区域,计算其顶点坐标在目标投影坐标系下的同名点坐标;
S205、为每一个所述投影变换区域,建立投影变换线性方程组;
S206、对于每一个所述投影变换区域,分别根据其顶点坐标以及相应的同名点坐标,计算对应的投影变换线性方程组的变换系数的值;
S207、将每一个所述投影变换区域的编码以及对应的变换系数、顶点坐标、同名点坐标进行存储;
S208、根据待投影要素的经度和纬度得到其对应的所述投影转换区域的编码,根据所述编码得到所述待投影要素所在的所述投影变换区域的顶点坐标、同名点坐标以及变换系数,得到所述待投影要素所在的投影变换区域的投影变换线性方程组;
S209、利用所述步骤S208得到的投影变换线性方程组计算所述待投影要素在目标投影坐标系下的坐标。
优选地,所述投影变换区域为矩形;对于每一个所述投影变换区域分别提取其四个顶点坐标:A(Ba,La)、B(Bb,Lb)、C(Bc,Lc)、D(Bd,Ld),其中A点和B点同一条纬线上,C点和D点在同一条纬线上,A点和D点在同一条经线上,B点和C点在同一条经线上;
对于每一个所述投影变换区域,分别计算其顶点A、B、D在目标投影坐标系下的同名点坐标A’(xa,ya),B’(xb,yb),D’(xd,yd)。
优选地,所述投影变换线性方程组为:
y=yaB(B-Ba)+θL(L-La)
其中(B,L),(x,y)为待投影要素分别在地理坐标系与目标投影坐标系下的坐标;θL,θB分别为所述变换系数;xa,xb分别为对应的所述投影变换区域的同名点A’的坐标值;Ba、La分别为对应的所述投影变换区域的顶点A的坐标值。
优选地,每一个所述投影变换区域的所述变换系数的求解方法均为:
对于一所述投影变换区域,将其点A、A’、B、B’、D、D’的坐标值带入对应的所述投影变换线性方程组,建立所述变换系数的计算公式:
θL=(yb-ya)(Bd-Ba)-(yd-ya)(Bb-Ba)/(Bd-Ba)(Lb-La)-(Bb-Ba)(Ld-La)
θB=(yd-ya)(Lb-La)-(yb-ya)(Ld-La)/(Bd-Ba)(Lb-La)-(Bb-Ba)(Ld-La)
将等值关系Ba=Bb,La=Ld代入所述变换系数的计算公式,可得到简化的公式:
θL=(yb-ya)/(Lb-La)
θB=(yd-ya)/(Bd-Ba)
根据点A、A’、B、B’、D、D’坐标值计算所述变换系数。
优选地,所述步骤S201中,以固定的经度跨度和固定纬度跨度为步长,并且以经线与纬线为边界划分所述投影变换区域。
优选地,所述步骤S202中,编码方式为:以经度和纬度最小的所述投影变换区域为编码起始点,编码规则为行号-列号,纬度范围相同的述投影变换区域为相同的行号,经度范围相同的所述投影变换区域为相同的列号。
优选地,所述编码的行号和列号由下面公式计算得到:
其中NB为所述编码的行号,NL为所述编码的列号,θ、γ分别为所述待投影要素的经度、纬度;Δb表示所述纬度跨度,Δl表示所述经度跨度。
(三)有益效果
本发明提供了地理坐标投影变化方法,所述方法将投影区域进行处理,得到投影变换区域,在很小的投影变换区域内将经线与纬线分别用直线表示,将球面计算转化为平面计算,可以提高地理坐标投影变换的计算速度,并且可以通过投影变换区域的大小实现误差控制,能够以较少的成本满足涉及大量数据变换并且对变换精度与变换速度同时有较高要求的实际应用需求;此外,本发明通过将多个投影变换区域的顶点坐标和变换系数等数据进行存储,在使用时进行查询的方式,在较大范围内实现本发明的方法的工程化应用。
附图说明
为了更清楚地说明本发明实施例或现有技术中的技术方案,下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍,显而易见地,下面描述中的附图仅仅是本发明的一些实施例,对于本领域普通技术人员来讲,在不付出创造性劳动的前提下,还可以根据这些附图获得其他的附图。
图1A为本发明的一种地理坐标投影变化方法的流程图;
图1B为本发明的一种地理坐标投影变化方法的投影变换区域的顶点分布示意图;
图2A为本发明的另一种地理坐标投影变化方法的流程图;
图2B为本发明的另一种地理坐标投影变化方法的投影变换区域划分方案示意图;
图2C为本发明的另一种地理坐标投影变化方法的投影变换区域编码规则示意图。
具体实施方式
下面结合附图和实施例对本发明作进一步详细描述。以下实施例用于说明本发明,但不能用来限制本发明的范围。
图1A为本发明的一种地理坐标投影变化方法的流程图;所述方法包括以下步骤:
S101、将地理坐标系下的投影区域规划为一个二维的投影变换区域,所述投影变换区域为能够覆盖所述投影区域的面积最小的几何图形;优选地为矩形,所述投影变换区域为能够覆盖所述投影区域的面积最小的矩形;
地球可以十分近似地模拟为一个旋转椭球体,这就是说,地球的各个经线分别是同一个南北对称轴、具有相同长半轴与短半轴、并在赤道平面上均匀立交分布的椭圆,而各个纬线则是半径不等的同心圆。由于经线圈十分巨大,而纬线圈即使在低纬度地区也非常巨大,这两种线圈的很小局部区域投影在地图平面上都可以视为直线。因此,在很小局部区域可以将经线与纬线分别用直线表示,利用经线和纬线划分投影变换区域或将不规则的、小的投影区域规划为矩形的投影变换区域,并将投影变换区域视为平面,将经纬度与平面坐标之间的复杂函数关系转化为简单的线性关系。但是这种线性转换的精度会随着纬度的增大以及变换区域范围的增大而降低,在具体实施中,需要根据投影变换的地理位置以及实际变换的精度要求设定坐标投影变换范围,例如在北纬20°-25°工作区域,设定坐标投影变换范围为3′*3′,可以保证变换误差在2米以内;
S102、提取所述变换区域的至少三个顶点的坐标,对于矩形的投影变换区域,提取其四个顶点坐标:A(Ba,La)、B(Bb,Lb)、C(Bc,Lc)、D(Bd,Ld),其中A点和B点同一条纬线上,C点和D点在同一条纬线上,A点和D点在同一条经线上,B点和C点在同一条经线上,如图1B所示;坐标值精确到小数点后7位,并且由于变换区域由经线和纬线划分,存在关系Ba=Bb,La=Ld
S103、计算顶点坐标在目标投影坐标系下的同名点坐标;对于矩形的投影变换区域,计算所述顶点A、B、D在目标投影坐标系下的同名点坐标A’(xa,ya),B’(xb,yb),D’(xd,yd);所述目标投影坐标系即为地理坐标执行投影变换后的平面坐标系;
所述同名点坐标的计算过程是经典地图学已解决的问题,通用模型包括在坐标基准不变情况下的解析法、数值法以及坐标基准变化情况下的相似变换模型等,理论上变换精度完全可以满足实际要求。因此,此步骤参照通用模型执行同名点坐标的高精度计算,并在接下来的算法流程中将同名点坐标的计算结果视为“真值”,也就是说同名点平面坐标计算精度是影响本发明的方法的计算精度的关键因素之一;
S104、设A(Ba,La)为原点,则经度与纬度相对于原点A的增量与平面坐标XY存在简单的线性对应关系,建立投影变换线性方程组:
y=yaB(B-Ba)+θL(L-La)
其中(B,L),(x,y)为待投影要素分别在地理坐标系与目标投影坐标系下的坐标;θL,θB分别为所述变换系数;其中θB为主变换系数,起着建立两者比例关系的重要作用,θL为辅助变换系数,调整经线方向增量、纬线方向增量;需要强调的是所述投影变换线性方程组的变换系数并不是固定的,其因原点A地理坐标、其同名点坐标的不同而不同,即因投影变换区域的位置以及大小的不同而不同;
S105、根据所述顶点A、B、D的顶点坐标以及相应的同名点坐标,计算所述投影变换线性方程组的变换系数的值:
将点A、A’、B、B’、D、D’的坐标值带入所述投影变换线性方程组,建立所述变换系数的计算公式:
θL=(yb-ya)(Bd-Ba)-(yd-ya)(Bb-Ba)/(Bd-Ba)(Lb-La)-(Bb-Ba)(Ld-La)
θB=(yd-ya)(Lb-La)-(yb-ya)-(yb-ya)/(Bd-Ba)(Lb-La)-(Bb-Ba)(Ld-La)
将等值关系Ba=Bb,La=Ld代入所述变换系数的计算公式,可得到简化的公式:
θL=(yb-ya)/(Lb-La)
θB=(yd-ya)/(Bd-Ba)
根据点A、A’、B、B’、D、D’坐标值,应用上面的公式计算所述变换系数;
S106、利用所述投影变换线性方程组计算所述投影变换区域内要素的目标投影坐标系下的坐标;以在投影变换区域内要素F进行快速投影变换为例,提取F点地理坐标(Bf,Lf)并带入步骤104建立的投影变换线性方程组,计算所述投影变换区域内对应的平面坐标XY。
本发明提出的方法可以在投影变换区域实现经纬度数据的快速投影变换,但是所述投影变换线性方程组计算过程必须限制变换区域范围以保障变换精度。
图2A为本发明的另一种地理坐标投影变化方法的流程图;本实施例以地理坐标系下北纬21°~24°,东经120°~123°的空间范围为工作区域对要素执行投影变换:
S201、将投影区域划分为两个或两个以上的二维投影变换区域,优选地为矩形,将整个投影区域在地理坐标系下视为平面,则经线与纬线为相互垂直的直线;在经线和纬线方向均以3′(3分)为步长,将工作区域划分为60*60个以经线与纬线为边界的投影变换区域,如图2B所示;
S202、为每个所述投影变换区域设置编码;如图2C所示,所述编码方式为:以经度和纬度最小的所述投影变换区域为编码起始点,编码规则为行号-列号,如点208所示,纬度范围相同的述投影变换区域为相同的列号,经度范围相同的所述投影变换区域为相同的行号;如点205所示,在206方向上对地理网格定义逐一递增的行号,在207方向上对地理网格定义逐一递增的列号;将所述待投影要素的经度和纬度转换为十进制分形式,
其中NB为对应的行号,NL为对应的列号,θ、γ分别为所述待投影要素的经度、纬度;Δb表示所述纬度跨度,Δl表示所述经度跨度。“7200”为东经120°用“分”表示的数字,“1260”为北纬21°用“分”表示的数字,“[]”为取整算符。“NB_NL”即为待投影要素所属投影变换区域的编码;
S203、对于每一个所述投影变换区域,分别提取其至少3个顶点坐标;对于矩形的投影变换区域,提取其四个顶点坐标:A(Ba,La)、B(Bb,Lb)、C(Bc,Lc)、D(Bd,Ld),其中A点和B点同一条纬线上,C点和D点在同一条纬线上,A点和D点在同一条经线上,B点和C点在同一条经线上;
S204对于每一个所述投影变换区域,计算其顶点坐标在目标投影坐标系下的同名点坐标;于矩形的投影变换区域,分别计算其顶点A、B、D在目标投影坐标系下的同名点坐标A’(xa,ya),B’(xb,yb),D’(xd,yd);
S205、为每一个所述投影变换区域,建立投影变换线性方程组;
S206、对于每一个所述投影变换区域,分别根据其顶点A、B、D的顶点坐标以及相应的同名点坐标,计算对应的投影变换线性方程组的变换系数的值;
上述步骤S203-S206与S102-S105中计算变换系数的方法相同,其中变换系数的计算过程不会受到相邻投影变换区域的影响;
S207、将所述投影变换区域的编码以及对应的变换系数、所述投影变换区域的顶点坐标、同名点坐标进行存储,存储到数据表中;表1为线性方程组参数表,通过编码可以快速检索到投影变换线性方程组的变换系数以及原点坐标,得到对应的投影变换线性方程组;
S208、根据待投影要素的经度和纬度将待投影要素的地理坐标转换为十进制分表示,得到NB、NL,查询得到其对应的所述投影转换区域的编码,根据所述编码得到所述待投影要素所在的所述投影变换区域的顶点坐标、同名点坐标以及变换系数,从而得到所述待投影要素所在的投影变换区域的投影变换线性方程组;在获取每个投影变换区域的变换系数以及原点坐标的基础上,可以针对每个投影变换区域建立变换系数已知投影变换线性方程组。但是由于针对不同投影变换区域建立的投影变换线性方程组并不统一,在待投影转换要素执行投影变换前需要确定所述要素所属的投影转换区域,以确定投影变换线性方程组的变换系数,这一计算过程也会影响到整个工程化投影变换措施的计算效率;
S209、利用所述步骤S208得到的投影变换线性方程组计算所述待投影要素在目标投影坐标系下的坐标。
表1
采用本发明的方法,可以提高地理坐标投影变换的计算速度,并且可以通过控制投影变换区域的大小实现误差控制,能够以较少的成本满足某些涉及大量数据变换并且对变换精度与变换速度同时有较高要求的实际应用需求;并且可以通过将较大范围的工作区域划分为相互独立的投影变换区域来实施本发明的方法的工程化应用。
以上实施方式仅用于说明本发明,而非对本发明的限制。尽管参照实施例对本发明进行了详细说明,本领域的普通技术人员应当理解,对本发明的技术方案进行各种组合、修改或者等同替换,都不脱离本发明技术方案的精神和范围,均应涵盖在本发明的权利要求范围当中。

Claims (9)

1.一种地理坐标投影变化方法,其特征在于,包括以下步骤:
S101、将地理坐标系下的投影区域规划为一个二维的投影变换区域,所述投影变换区域为能够覆盖所述投影区域的面积最小的几何图形;
S102、提取所述投影变换区域的至少三个顶点的坐标;
S103、计算步骤S102的顶点坐标在目标投影坐标系下的同名点坐标;
S104、建立投影变换线性方程组,所述投影变换线性方程组为:
y=yaB(B-Ba)+θL(L-La)
其中(B,L),(x,y)为待投影要素分别在地理坐标系与目标投影坐标系下的坐标;θL,θB分别为所述变换系数;
S105、根据所述顶点坐标以及相应的同名点坐标,计算所述投影变换线性方程组的变换系数的值;
S106、利用所述投影变换线性方程组计算所述投影变换区域内要素在目标投影坐标系下的坐标。
2.根据权利要求1所述的方法,其特征在于,所述投影变换区域为矩形,提取其四个顶点的顶点坐标:A(Ba,La)、B(Bb,Lb)、C(Bc,Lc)、D(Bd,Ld),其中A点和B点同一条纬线上,C点和D点在同一条纬线上,A点和D点在同一条经线上,B点和C点在同一条经线上;计算所述顶点A、B、D在目标投影坐标系下的同名点坐标A’(xa,ya),B’(xb,yb),D’(xd,yd)。
3.根据权利要求2所述的方法,其特征在于,所述变换系数的求解方法具体为:
将点A、A’、B、B’、D、D’的坐标值带入所述投影变换线性方程组,建立所述变换系数的计算公式:
θL=(yb-ya)(Bd-Ba)-(yd-ya)(Bb-Ba)/(Bd-Ba)(Lb-La)-(Bb-Ba)(Ld-La)
θB=(yd-ya)(Lb-La)-(yb-ya)(Ld-La)/(Bd-Ba)(Lb-La)-(Bb-Ba)(Ld-La)
将等值关系Ba=Bb,La=Ld代入所述变换系数的计算公式,可得到简化的公式:
θL=(yb-ya)/(Lb-La)
θB=(yd-ya)/(Bd-Ba)
根据点A、A’、B、B’、D、D’坐标值计算所述变换系数。
4.一种地理坐标投影变化方法,其特征在于,包括以下步骤:
S201、将投影区域划分为两个或两个以上的二维投影变换区域;
S202、为每个所述投影变换区域设置编码;
S203、对于每一个所述投影变换区域,提取其至少三个顶点的坐标;
S204、对于每一个所述投影变换区域,计算其顶点坐标在目标投影坐标系下的同名点坐标;
S205、为每一个所述投影变换区域,建立投影变换线性方程组,所述投影变换线性方程组为:
y=yaB(B-Ba)+θL(L-La)
其中(B,L),(x,y)为待投影要素分别在地理坐标系与目标投影坐标系下的坐标;θL,θB分别为所述变换系数;xa,xb分别为对应的所述投影变换区域的同名点A’的坐标值;Ba、La分别为对应的所述投影变换区域的顶点A的坐标值;
S206、对于每一个所述投影变换区域,分别根据其顶点坐标以及相应的同名点坐标,计算对应的投影变换线性方程组的变换系数的值;
S207、将每一个所述投影变换区域的编码以及对应的变换系数、顶点坐标、同名点坐标进行存储;
S208、根据待投影要素的经度和纬度得到其对应的所述投影转换区域的编码,根据所述编码得到所述待投影要素所在的所述投影变换区域的顶点坐标、同名点坐标以及变换系数,得到所述待投影要素所在的投影变换区域的投影变换线性方程组;
S209、利用所述步骤S208得到的投影变换线性方程组计算所述待投影要素在目标投影坐标系下的坐标。
5.根据权利要求4所述的方法,其特征在于,所述投影变换区域为矩形;对于每一个所述投影变换区域分别提取其四个顶点坐标:A(Ba,La)、B(Bb,Lb)、C(Bc,Lc)、D(Bd,Ld),其中A点和B点同一条纬线上,C点和D点在同一条纬线上,A点和D点在同一条经线上,B点和C点在同一条经线上;
对于每一个所述投影变换区域,分别计算其顶点A、B、D在目标投影坐标系下的同名点坐标A’(xa,ya),B’(xb,yb),D’(xd,yd)。
6.根据权利要求4所述的方法,其特征在于,每一个所述投影变换区域的所述变换系数的求解方法均为:
对于一所述投影变换区域,将其点A、A’、B、B’、D、D’的坐标值带入对应的所述投影变换线性方程组,建立所述变换系数的计算公式:
θL=(yb-ya)(Bd-Ba)-(yd-ya)(Bb-Ba)/(Bd-Ba)(Lb-La)-(Bb-Ba)(Ld-La)
θB=(yd-ya)(Lb-La)-(yb-ya)(Ld-La)/(Bd-Ba)(Lb-La)-(Bb-Ba)(Ld-La)
将等值关系Ba=Bb,La=Ld代入所述变换系数的计算公式,可得到简化的公式:
θL=(yb-ya)/(Lb-La)
θB=(yd-ya)/(Bd-Ba)
根据点A、A’、B、B’、D、D’坐标值计算所述变换系数。
7.根据权利要求6所述的方法,其特征在于,所述步骤S201中,以固定的经度跨度和固定纬度跨度为步长,并且以经线与纬线为边界划分所述投影变换区域。
8.根据权利要求7所述的方法,其特征在于,所述步骤S202中,编码方式为:以经度和纬度最小的所述投影变换区域为编码起始点,编码规则为行号-列号,纬度范围相同的述投影变换区域为相同的行号,经度范围相同的所述投影变换区域为相同的列号。
9.根据权利要求8所述的方法,其特征在于,所述编码的行号和列号由下面公式计算得到:
其中NB为所述编码的行号,NL为所述编码的列号,θ、γ分别为所述待投影要素的经度、纬度;Δb表示所述纬度跨度,Δl表示所述经度跨度。
CN201410750917.3A 2014-12-09 2014-12-09 地理坐标投影变化方法 Expired - Fee Related CN104537606B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201410750917.3A CN104537606B (zh) 2014-12-09 2014-12-09 地理坐标投影变化方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201410750917.3A CN104537606B (zh) 2014-12-09 2014-12-09 地理坐标投影变化方法

Publications (2)

Publication Number Publication Date
CN104537606A CN104537606A (zh) 2015-04-22
CN104537606B true CN104537606B (zh) 2018-02-02

Family

ID=52853125

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201410750917.3A Expired - Fee Related CN104537606B (zh) 2014-12-09 2014-12-09 地理坐标投影变化方法

Country Status (1)

Country Link
CN (1) CN104537606B (zh)

Families Citing this family (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105243268B (zh) * 2015-09-18 2018-06-22 网易(杭州)网络有限公司 一种游戏地图定位方法、装置及用户终端
CN108287915A (zh) * 2018-02-11 2018-07-17 浙江科澜信息技术有限公司 地理信息系统中测绘成果的坐标转换方法及系统
CN108681986B (zh) * 2018-04-28 2020-03-31 吉林大学 一种加快正弦投影格式遥感产品投影变换的方法及系统
CN109191579B (zh) * 2018-07-02 2019-06-28 广东省地质测绘院 坐标转换的方法、装置、计算机设备和存储介质
CN109300391B (zh) * 2018-10-16 2020-11-03 中国人民解放军海军工程大学 一种地图设计方法和装置
CN109974716B (zh) * 2019-02-27 2021-04-06 上海晶赞融宣科技有限公司 地理位置匹配方法及装置、存储介质、服务器
CN111198927B (zh) * 2020-01-02 2020-11-17 平安科技(深圳)有限公司 地理位置数据的匹配方法、装置、电子设备及存储介质
CN112102430B (zh) * 2020-08-18 2023-11-21 国家海洋信息中心 单要素变化曲线至地理地图的投影计算方法

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101833896A (zh) * 2010-04-23 2010-09-15 西安电子科技大学 基于增强现实的地理信息指引方法与系统
CN102520406A (zh) * 2011-12-20 2012-06-27 中国测绘科学研究院 基于影像坐标映射的sar立体模型构建与量测方法

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
KR20040061244A (ko) * 2002-12-30 2004-07-07 삼성전자주식회사 디-인터레이싱 방법 및 그 장치

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101833896A (zh) * 2010-04-23 2010-09-15 西安电子科技大学 基于增强现实的地理信息指引方法与系统
CN102520406A (zh) * 2011-12-20 2012-06-27 中国测绘科学研究院 基于影像坐标映射的sar立体模型构建与量测方法

Also Published As

Publication number Publication date
CN104537606A (zh) 2015-04-22

Similar Documents

Publication Publication Date Title
CN104537606B (zh) 地理坐标投影变化方法
CN103389080B (zh) 一种基于地理信息应用城市独立坐标系参数获取方法
CN107193888B (zh) 一种面向车道级导航定位的城市道路路网模型
CN105160197A (zh) 一种综合性地理空间数据坐标转换方法及系统
CN105203104A (zh) 一种适用于高精度惯导系统的重力场建模方法
CN109801219A (zh) 面向在线地图叠加的gis数据校正方法及装置
WO2019062166A1 (zh) 跨平台月基对地观测影像自动几何校正方法
Kennedy Mathematical geography
CN104463970A (zh) 一种基于遥感影像确定城市三维重心的方法及其应用
CN110866015A (zh) 一种基于局部网格的动目标活动范围记录方法
CN103530904A (zh) 一种基于克里金方法的水下地形数字高程建立方法
CN107462220B (zh) 面向月基对地观测影像的投影极坐标几何表达方法
CN104700710A (zh) 一种房产测绘仿真地图
CN114580310A (zh) 一种基于palm实现wrf模拟风场降尺度处理的方法
CN114881466A (zh) 一种基于多源数据的人口空间分区拟合方法
CN104932858A (zh) 一种精确匹配实时漫游的背景地图叠加显示方法
CN114119903A (zh) 一种基于实景三维城市的动态交通模拟方法
CN204988280U (zh) 简便高效平面任意线路工程测设装置
CN114187399A (zh) 一种地外星体表面光照仿真影像生成方法和系统
CN114387419A (zh) 一种基于多维空间剖分的三维地理实体编码方法及装置
Moscovici et al. Mathematical cartography based on georeferencing maps
Setiawan et al. Calculation Of Central Java Province Region Area Using Shoelace Formula Based On The GADM Database
CN107220314A (zh) 一种构建数字高程模型的方法
CN111710016B (zh) 一种晨昏线图层获取方法及装置
Akin et al. “Yu laid out the lands:” georeferencing the Chinese Yujitu [Map of the Tracks of Yu] of 1136

Legal Events

Date Code Title Description
PB01 Publication
C10 Entry into substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant
CF01 Termination of patent right due to non-payment of annual fee

Granted publication date: 20180202

CF01 Termination of patent right due to non-payment of annual fee