CN110390620B - 一种基于曲线傅里叶变换的等高线地图数据盲水印方法 - Google Patents

Abstract

本发明公开了一种基于曲线傅里叶变换的等高线地图数据盲水印方法。传统水印算法还很难直接应用在等高线地图数据上。本发明方法包括水印信息嵌入和水印信息提取。水印信息嵌入是把用户选择的水印信息嵌入到原始等高线地图中得到一个新的等高线地图，具体包括调制水印信息、对等高线地图数据进行处理和嵌入水印信息。水印信息提取是从待检测的等高线地图中提取水印信息内容，具体包括对待检测的等高线地图数据进行处理、提取水印信息。本发明具有较好的鲁棒性,并且是一种盲检方法,具有较好的实用价值。

Description

一种基于曲线傅里叶变换的等高线地图数据盲水印方法

技术领域

本发明属于多媒体信息安全技术领域，涉及一种针对地理空间数据中等高线地图数据的盲水印方法。

背景技术

数字水印技术是将代表著作权人身份的特定信息，按照某种方式嵌入被保护的信息中,在产生版权纠纷时,通过相应的算法提取出该数字水印，从而验证版权的归属，确保著作权人的合法利益，避免非法盗版的威胁。被保护的信息可以是任何一种数字媒体,如图像、声音、视频或一般性的电子文档等。数字水印技术运用于数据版权保护和数据真伪检测应用中，是一种可在开放的网络环境下实现信息隐藏与跟踪的新型技术。

等高线地图数据水印算法不同于传统的基于图像水印算法，这主要是由于等高线数据不具有图像数据应用的普遍性等原因,同时还在于等高线地图数据在数据结构、存储形式、表现方式、应用环境、使用要求以及可能的攻击行为等方面与图像等多媒体数据都不一样,一般的数字水印算法还很难直接应用在等高线地图数据上，这使得目前针对等高线地图的数字水印技术应用还非常匮乏。

由于等高线地图是地理空间数据的一个重要组成部分，它是描述关于人类赖以生存的地球的重要信息，是国家基础性、公益性设施，是国家基础设施建设和地球科学研究的支撑性成果，是国家经济、国防建设中不可缺少的资源，对社会持续发展起着重要作用，因此对于等高线地图数据的水印技术也具有非常好的使用前景。由此可见，数字水印技术在等高线地图数据上的应用具有更加潜在的用途，而且具有非常好的经济社会效益。

发明内容

本发明针对现有技术的不足，提供了一种基于曲线傅里叶变换的等高线地图数据盲水印方法。

本发明方法的包括水印信息嵌入和水印信息提取。

所述的水印信息嵌入是把用户选择的水印信息嵌入到原始等高线地图中得到一个新的等高线地图，输入为一个等高线地图M和一个包含11个0～9之间整数的序列w，输出为包含水印信息的等高线地图M′，具体步骤为：

步骤(1):曲线傅里叶变换计算

设有一顶点序列D₁，D₂，...，D_i，其对应的坐标序列为(x₁,y₁)，(x₂,y₂)，...，(x_i,y_i)，对于顶点序列D₁，D₂，...，D_i，使用曲线傅里叶变换计算得到一组曲线傅里叶变换系数a_n，b_n，c_n和d_n。

步骤(2):调制水印信息

嵌入的水印信息是一个包含11个整数的序列，每个整数取值在1到9之间，记这个整数序列为w＝(w₁，w₂，...，w_m，...，w₁₁)。

画一圆R，其半径为1，在圆R中分别画角度为0度、30度、60度，90度、120度，150度，180度、210度、240度，270度、300度、330度的半径线，这12条半径线与R的交点分别记为A₀、A₁、A₂、A₃、A₄、A₅、A₆、A₇、A₈、A₉、A₁₀及A₁₁。其中，A₀是0度半径线与R的交点，A₁是30度半径线与的R的交点，A₃是60度半径线与R的交点，其他以此类推。

把A_k交点在圆R上的坐标分别记为(x_k,y_k),k＝0,1,，...，11那么有：

x_k＝cos(k×π/6)，y_k＝sin(k×π/6)

将w＝(w₁，w₂，...，w_m，...，w₁₁)的11个整数值按从前到后的顺序，分别与A₁、A₂、A₃、A₄、A₅、A₆、A₇、A₈、A₉、A₁₀、A₁₁一一对应，然后把A₁、A₂、A₃、A₄、A₅、A₆、A₇、A₈、A₉、A₁₀、A₁₁的横坐标和纵坐标分别除以w₁，w₂，w₃，w₄，w₅，w₆，w₇，w₈，w₉，w₁₀，w₁₁的整数值。

计算A₀和点A₁之间线段长度h₀；计算A₁和点A₂之间线段长度h₁，计算A₂和A₃之间的线段长度h₂、A₃和A₄之间的线段长度h₃、A₄和A₅之间的线段长度h₄、A₅和A₆之间的线段长度h₅、A₆和A₇之间的线段长度h₆、A₇和A₈之间的线段长度h₇、A₈和A₉之间的线段长度h₈、A₉和A₁₀之间的线段长度h₉、A₁₀和A₁₁之间的线段长度h₁₀及A₁₁和A₀之间的线段长度h₁₁。

对A₀、A₁、A₂、A₃、A₄、A₅、A₆、A₇、A₈、A₉、A₁₀及A₁₁使用曲线傅里叶变换得到一组曲线傅里叶变换系数，得到系数Wa_n、Wb_n、Wc_n和Wd_n。这样，根据水印信息w＝(w₁，w₂，...，w_m，...，w₁₁)，利用曲线傅里叶变换得到一组曲线傅里叶变换系数Wa_n，Wb_n，Wc_n，Wd_n，n＝1，2，...，199，200。

步骤(3)：嵌入水印信息

步骤(a)：

对等高线地图的所有闭合曲线分别计算其周长，对所有闭合曲线按周长进行排序，周长越大，排序越靠前。

对等高线地图的所有闭合曲线排好序后，在排好序的闭合曲线中选择最靠前的五条闭合曲线作为嵌入数字水印的对象。将排序第一的闭合曲线记为L₁,排序第二的闭合曲线记为L₂，排序第三，四，五的闭合曲线分别记为L₃，L₄和L₅。

将L₁的周长记为TL₁，L₂的周长记为TL₂，L₃的周长记为TL₃，L₄的周长记为TL₄，L₅的周长记为TL₅。

步骤(b)：

对L₁，从0度水平开始按逆时针顺序扫描闭合曲线所记录包含的顶点，获得一个顶点序列。对于该顶点序列，使用曲线傅里叶变换得到一组曲线傅里叶变换系数，分别记为L₁a_n、L₁b_n、L₁c_n和L₁d_n，其中n＝1，2，...，999，1000。

在使用曲线傅里叶变换计算L₁顶点序列的曲线傅里叶系数时，同时保存所有的t_j，也就是保存L₁第一个顶点沿着曲线按逆时针方向到其他每一个顶点的曲线距离。

采用与如上L₁同样的计算方法，对于L₂计算得到一组曲线傅里叶变换系数L₂a_n，L₂b_n，L₂c_n，L₂d_n，n＝1，2，...，999，1000和L₂第一个顶点沿着曲线按逆时针方向到其他每一个顶点的曲线距离。

采用与如上L₁同样的方法，对于L₃计算得到一组曲线傅里叶变换系数L₃a_n，L₃b_n，L₃c_n，L₃d_n，n＝1，2，...，999，1000和L₃第一个顶点沿着曲线按逆时针方向到其他每一个顶点的曲线距离。

采用与如上L₁同样的方法，对于L₄计算得到一组曲线傅里叶变换系数L₄a_n，L₄b_n，L₄c_n，L₄d_n，n＝1，2，...，999，1000和L₄第一个顶点沿着曲线按逆时针方向到其他每一个顶点的曲线距离。

采用与如上L₁同样的方法，对于L₅计算得到一组曲线傅里叶变换系数L₅a_n，L₅b_n，L₅c_n，L₅d_n，n＝1，2，...，999，1000和L₅第一个顶点沿着曲线按逆时针方向到其他每一个顶点的曲线距离。

步骤(c):

计算

将Wa_n，Wb_n，Wc_n，Wd_n，n＝1，2，...，50的这些系数都分别乘以

计算

将Wa_n，Wb_n，Wc_n，Wd_n，n＝51，52，...，100的这些系数都分别乘以

计算

将Wa_n，Wb_n，Wc_n，Wd_n，n＝101，102，...，150的这些系数都分别乘以

计算

将Wa_n，Wb_n，Wc_n，Wd_n，n＝151，152，...，200的这些系数都分别乘以

计算

将h₀，h₁，h₂，h₃，h₄，h₅，h₆，h₇，h₈，h₉，h₁₀，h₁₁都分别乘以

其中[·]表示取整运算。

步骤(d):

将Wa_n，n＝1，2，...，50这些系数按n的顺序放在L₁a₁₀₀₀系数后面，得到L₁a₁，L₁a₂，...，L₁a₁₀₀₀，Wa₁，Wa₂，...，Wa₅₀序列。将Wa_n，n＝51，52，...，100这些系数按n的顺序放在L₂a₁₀₀₀系数后面；将Wa_n，n＝101，102，...，150这些系数按n的顺序放在L₃a₁₀₀₀系数后面；将Wa_n，n＝151，152，...，200这些系数按n的顺序放在L₄a₁₀₀₀系数后面。

将Wb_n，n＝1，2，...，50这些系数按n的顺序放在L₁b₁₀₀₀系数后面，得到L₁b₁，L₁b₂，...，L₁b₁₀₀₀，Wb₁，Wb₂，...，Wb₅₀序列；将Wb_n，n＝51，52，...，100这些系数按n的顺序放在L₂b₁₀₀₀系数后面；将Wb_n，n＝101，102，...，150这些系数按n的顺序放在L₃b₁₀₀₀系数后面；将Wb_n，n＝151，152，...，200这些系数按n的顺序放在L₄b₁₀₀₀系数后面。

同样，使用上述Wa_n和Wb_n分片放置的方式，将Wc_n的系数分片按顺序放到L₁c₁₀₀₀，L₂c₁₀₀₀，L₃c₁₀₀₀，L₄c₁₀₀₀的后面。将Wd_n的系数分片按顺序放到L₁d₁₀₀₀，L₂d₁₀₀₀，L₃d₁₀₀₀，L₄d₁₀₀₀的后面。

将h₀，h₁，h₂，h₃，h₄，h₅，h₆，h₇，h₈，h₉，h₁₀，h₁₁按顺序分别放在L₅a₁₀₀₀，L₅b₁₀₀₀，L₅c₁₀₀₀和L₅d₁₀₀₀的后面，分别得到L₅a₁，L₅a₂，...，L₅a₁₀₀₀，h₀，h₁，...，h₁₀，h₁₁序列；L₅b₁，L₅b₂，...，L₅b₁₀₀₀，h₀，h₁，...，h₁₀，h₁₁序列；L₅c₁，L₅c₂，...，L₅c₁₀₀₀，h₀，h₁，...，h₁₀，h₁₁序列；L₅d₁，L₅d₂，...，L₅d₄₉，L₅d₁₀₀₀，h₀，h₁，...，h₁₀，h₁₁序列。

步骤(e):

对于L₁a₁，L₁a₂，...，L₁a₁₀₀₀，Wa₁，Wa₂，...，Wa₅₀；L₁b₁，L₁b₂，...，L₁b₁₀₀₀，Wb₁，Wb₂，...，Wb₅₀；L₁c₁，L₁c₂，...，L₁c₅₀，Wc₁，Wc₂，...，Wc₁₀₀₀；L₁d₁，L₁d₂，...，L₁d₁₀₀₀，Wd₁，Wd₂，...，Wd₅₀，采用逆曲线傅里叶变换得到顶点坐标。

由顶点坐标得到一个顶点序列(x(n),y(n))，n＝0，1，2，...，N-1。将这个顶点序列依次用线段连接得到一个闭合曲线WL₁。把这个闭合曲线WL₁替代L₁保存到等高线地图中。

对于L₂a₁，L₂a₂，...，L₂a₁₀₀₀，Wa₅₁，Wa₅₂，...，Wa₁₀₀；L₂b₁，L₂b₂，...，L₂b₁₀₀₀，Wb₅₁，Wb₅₂，...，Wb₁₀₀；L₂c₁，L₂c₂，...，L₂c₁₀₀₀，Wc₅₁，Wc₅₂，...，Wc₁₀₀；L₂d₁，L₂d₂，...，L₂d₁₀₀₀，Wd₅₀，Wd₅₁，...，Wd₁₀₀，采用与WL₁同样的方法，利用逆曲线傅里叶变换得到一个顶点序列，将这个顶点序列依次连接得到WL₂，将WL₂替代L₂保存到等高线地图中。

对于L₃a₁，L₃a₂，...，L₃a₁₀₀₀，Wa₁₀₁，Wa₁₀₂，...，Wa₁₅₀；L₃b₁，L₃b₂，...，L₃b₁₀₀₀，Wb₁₀₁，Wb₁₀₂，...，Wb₁₅₀；L₃c₁，L₃c₂，...，L₃c₁₀₀₀，Wc₁₀₁，Wc₁₀₂，...，Wc₁₅₀；L₃d₁，L₃d₂，...，L₃d₁₀₀₀，Wd₁₀₁，Wd₁₀₂，...，Wd₁₅₀，采用与WL₁同样的方法，利用逆曲线傅里叶变换得到一个顶点序列，将这个顶点序列依次连接得到WL₃，将WL₃替代L₃保存到等高线地图中。

对于L₄a₁，L₄a₂，...，L₄a₁₀₀₀，Wa₁₅₁，Wa₁₅₂，...，Wa₂₀₀；L₄b₁，L₄b₂，...，L₄b₅₀，Wb₁₅₁，Wb₁₅₂，...，Wb₂₀₀；L₄c₁，L₄c₂，...，L₄c₁₀₀₀，Wc₁₅₁，Wc₁₅₂，...，Wc₂₀₀；L₄d₁，L₄d₂，...，L₄d₁₀₀₀，Wd₁₅₁，Wd₁₅₂，...，Wd₂₀₀，采用与WL₁同样的方法，利用逆曲线傅里叶变换得到一个顶点序列，将这个顶点序列依次连接得到WL₄，将WL₄替代L₄保存到等高线地图中。

对于，L₅a₁，L₅a₂，...，L₅a₁₀₀₀，h₀，h₁，...，h₁₀，h₁₁；L₅b₁，L₅b₂，...，L₅b₁₀₀₀，h₀，h₁，...，h₁₀，h₁₁；L₅c₁，L₅c₂，...，L₅c₁₀₀₀，h₀，h₁，...，h₁₀，h₁₁；L₅d₁，L₅d₂，...，L₅d₁₀₀₀，h₀，h₁，...，h₁₀，h₁₁，采用与WL₁同样的方法，利用逆曲线傅里叶变换得到一个顶点序列，将这个顶点序列依次连接可以得到WL₅，将WL₅替代L₅保存到等高线地图中。

这样，经过步骤(a)至步骤(e)，将水印信息w＝(w₁，w₂，...，w_m，...，w₁₁)嵌入到了等高线地图中。

所述的水印信息的提取是从待检测的等高线地图M′中提取水印信息w的内容，输入为包含水印信息的等高线地图M′，输出为等高线地图里包含的水印信息w，具体包括以下步骤：

步骤(1)：

对等高线地图的所有闭合曲线分别计算其周长，对所有闭合曲线按周长进行排序，周长越大，排序越靠前。

对等高线地图的所有闭合曲线排好序后，在排好序的闭合曲线中选择最靠前的五条闭合曲线作为潜入数字水印的对象。将排序第一的闭合曲线记为L₁,排序第二的闭合曲线记为L₂，排序第三，四，五的闭合曲线分别记为L₃，L₄和L₅。将L₁的周长记为TL₁，L₂的周长记为TL₂，L₃的周长记为TL₃，L₄的周长记为TL₄，L₅的周长记为TL₅。

步骤(2)：

步骤(a):

对L₁，从0度水平开始按逆时针顺序扫描闭合曲线所记录包含的顶点，获得一个顶点序列。对于该顶点序列，使用曲线傅里叶变换得到一组曲线傅里叶变换系数，记为L₁a_n、L₁b_n、L₁c_n和L₁d_n，n＝1，2，...，1050。

把L₁a₁₀₀₁赋值给Wa₁，把L₁a₁₀₀₂赋值给Wa₂，以此类推，这样得到Wa₁，Wa₂，...，Wa₅₀。把L₁b₁₀₀₁赋值给Wb₁，把L₁b₁₀₀₂赋值给Wb₂，以此类推，这样得到Wb₁，Wb₂，...，Wb₅₀。用同样的方法，利用L₁c_n和L₁d_n分别得到Wc₁，Wc₂，...，Wc₅₀和Wd₁，Wd₂，...，Wd₅₀。

步骤(b):

按如上所述步骤(a)的方法，对于L₂计算得到Wa_n，Wb_n，Wc_n，Wd_n，n＝51，52，...，100；对于L₃计算得到Wa_n，Wb_n，Wc_n，Wd_n，n＝101，102，...，150；对于L₄计算得到Wa_n，Wb_n，Wc_n，Wd_n，n＝151，152，...，200。

步骤(c):

对于L₅,采用步骤(a)L₁的方法得到顶点序列，然后对于该顶点序列，使用曲线傅里叶变换公式计算得到一组曲线傅里叶变换系数，n＝1，2，...，1012。最后把L₁a₁₀₀₁赋值给h₀，把L₁a₁₀₀₂赋值给h₁，以此类推，这样得到h₀，h₁，...，h₁₀，h₁₁。

步骤(d):

计算

将Wa_n，Wb_n，Wc_n，Wd_n，n＝1，2，...，50的这些系数都分别除以

计算

将Wa_n，Wb_n，Wc_n，Wd_n，n＝51，52，...，100，这些系数都分别除以

计算

将Wa_n，Wb_n，Wc_n，Wd_n，n＝101，102，...，150，这些系数都分别除以

计算

将Wa_n，Wb_n，Wc_n，Wd_n，n＝151，152，...，200，这些系数都分别除以

计算

将h₀，h₁，...，h₁₀，h₁₁这些系数都分别除以

步骤(3)：

按如下公式计算顶点坐标：

其中，n＝0，1，2，...，11。

计算

n＝1，2，...，11，得到一个实数序列。对实数序列L(n)，n＝1，2，...，11的每个实数分别进行四舍五入计算，然后得到一个Z(n)整数序列。这个整数序列就是提取出来的数字水印信息。

本发明的有益效果：本发明方法具有较好的鲁棒性,并且是一种盲检方法,具有较好的实用价值。

具体实施方式

一种基于曲线傅里叶变换的等高线地图数据盲水印方法，该方法包括水印信息嵌入和水印信息提取两部分。

(一)、所述的水印信息嵌入是把用户选择的水印信息嵌入到原始等高线地图中得到一个新的等高线地图，输入为一个等高线地图M和一个包含11个0～9之间整数的序列w，输出为包含水印信息的等高线地图M′，具体步骤包括：

步骤(1):曲线傅里叶变换计算

设有一顶点序列D₁，D₂，...，D_i，其对应的坐标序列为(x₁,y₁)，(x₂,y₂)，...，(x_i,y_i)，即D₁对应了(x₁,y₁)，D₂对应了(x₂,y₂)，后面以此类推。

对于顶点序列D₁，D₂，...，D_i，使用公式(1)到公式(4)所示的曲线傅里叶变换计算得到一组曲线傅里叶变换系数a_n，b_n，c_n，d_n。

其中，Δx_j＝x_j-x_j-1，Δy_j＝y_j-y_j-1，

T＝t_i，

n＝1,2,3...。实际上，t_j是第一个顶点沿着曲线到其他每个顶点的曲线长度。

步骤(2):调制水印信息

嵌入的水印信息w是一个包含11个整数的整数序列，每个整数取值在1到9之间，记这个整数序列为w＝(w₁，w₂，...，w_m，...，w₁₁)。

画一圆R，其半径为1，在圆R中分别画角度为0度、30度、60度，90度、120度，150度，180度、210度、240度，270度、300度、330度的半径线，这12条半径线与R的交点分别记为A₀、A₁、A₂、A₃、A₄、A₅、A₆、A₇、A₈、A₉、A₁₀及A₁₁。其中，A₀是0度半径线与R的交点，A₁是30度半径线与的R的交点，A₃是60度半径线与R的交点，其他以此类推。

把A_k交点在圆R上的坐标分别记为(x_k,y_k),k＝0,1,，...，11那么有：

x_k＝cos(k×π/6)，y_k＝sin(k×π/6),k＝0,1，...，11。

将w＝(w₁，w₂，...，w_m，...，w₁₁)的11个整数值按从前到后的顺序，分别与A₁、A₂、A₃、A₄、A₅、A₆、A₇、A₈、A₉、A₁₀、A₁₁一一对应，然后把A₁、A₂、A₃、A₄、A₅、A₆、A₇、A₈、A₉、A₁₀、A₁₁的横坐标和纵坐标分别除以和w₁，w₂，w₃，w₄，w₅，w₆，w₇，w₈，w₉，w₁₀，w₁₁的整数值。例如A₁对应了w₁，那么A₁的坐标调整为x₁＝cos(1×π/6)/w₁,y₁＝sin(1×π/6)/w₁。A₂对应了w₂，那么A₂的坐标值就调整为x₂＝cos(2×π/6)/w₂，y₂＝sin(2×π/6)/w₂，其他以此类推。

计算A₀和点A₁之间线段长度h₀，计算方法为

Δx₁＝x₁-x₀，y₁＝y₁-y₀。计算A₁和点A₂之间线段长度h₁，计算方法为

Δx₂＝x₂-x₁，y₂＝y₂-y₁。依次类推，计算A₂和A₃之间的线段长度h₂、A₃和A₄之间的线段长度h₃、A₄和A₅之间的线段长度h₄、A₅和A₆之间的线段长度h₅、A₆和A₇之间的线段长度h₆、A₇和A₈之间的线段长度h₇、A₈和A₉之间的线段长度h₈、A₉和A₁₀之间的线段长度h₉、A₁₀和A₁₁之间的线段长度h₁₀、及A₁₁和A₀之间的线段长度h₁₁。

对A₀、A₁、A₂、A₃、A₄、A₅、A₆、A₇、A₈、A₉、A₁₀及A₁₁使用公式(1)-公式(4)所示的曲线傅里叶变换公式计算，得到一组曲线傅里叶变换系数，其中n＝1，2，...，199，200。这里把公式(1)计算得到的系数记为Wa_n，把公式(2)计算得到的系数记为Wb_n，把公式(3)计算得到的系数记为Wc_n，把公式(4)计算得到的系数记为Wd_n。这样，根据水印信息w＝(w₁，w₂，...，w_m，...，w₁₁)，利用曲线傅里叶变换可以得到一组曲线傅里叶变换系数Wa_n，Wb_n，Wc_n，Wd_n，n＝1，2，...，199，200。

步骤(3)：对等高线地图M嵌入水印信息

水印信息的嵌入可以分为如下子步骤：

步骤(a)

对等高线地图M的所有闭合曲线分别计算其周长，对所有闭合曲线按周长进行排序，周长越大，排序越靠前。如果有曲线周长相同，那么对周长相等的闭合曲线再按顶点数量进行排序，顶点数量越大，排序越前。

对等高线地图M的所有闭合曲线排好序后，在排好序的闭合曲线中选择最靠前的五条闭合曲线作为嵌入数字水印的对象。这里将排序第一的闭合曲线记为L₁,排序第二的闭合曲线记为L₂，后面以此类推，排序第三，四，五的闭合曲线分别记为L₃，L₄和L₅。

将L₁的周长记为TL₁，L₂的周长记为TL₂，L₃的周长记为TL₃，L₄的周长记为TL₄，L₅的周长记为TL₅。

步骤(b):

对L₁，从0度水平开始按逆时针顺序扫描闭合曲线所记录包含的顶点，可以获得一个顶点序列。对于该顶点序列，使用公式(1)-公式(4)所示的曲线傅里叶变换公式计算得到一组曲线傅里叶变换系数，其中n＝1,2,...,999,1000。把公式(1)计算得到的系数记为L₁a_n，把公式(2)计算得到的系数记为L₁b_n，把公式(3)计算得到的系数记为L₁c_n，把公式(4)计算得到的系数记为L₁d_n。

在使用曲线傅里叶变换计算L₁顶点序列的曲线傅里叶系数时，同时保存所有的t_j，也就是保存L₁第一个顶点沿着曲线按逆时针方向到其他每一个顶点的曲线距离。

采用与如上L₁同样的计算方法，对于L₂计算得到一组曲线傅里叶变换系数L₂a_n，L₂b_n，L₂c_n，L₂d_n，n＝1,2,...,999,1000和L₂第一个顶点沿着曲线按逆时针方向到其他每一个顶点的曲线距离。

采用与如上L₁同样的方法，对于L₃计算得到一组曲线傅里叶变换系数L₃a_n，L₃b_n，L₃c_n，L₃d_n，n＝1,2,...,999,1000和L₃第一个顶点沿着曲线按逆时针方向到其他每一个顶点的曲线距离。

采用与如上L₁同样的方法，对于L₄计算得到一组曲线傅里叶变换系数L₄a_n，L₄b_n，L₄c_n，L₄d_n，n＝1,2,...,999,1000和L₄第一个顶点沿着曲线按逆时针方向到其他每一个顶点的曲线距离。

采用与如上L₁同样的方法，对于L₅计算得到一组曲线傅里叶变换系数L₅a_n，L₅b_n，L₅c_n，L₅d_n，n＝1,2,...,999,1000和L₅第一个顶点沿着曲线按逆时针方向到其他每一个顶点的曲线距离。

步骤(c):

计算

其中[·]表示取整运算。将Wa_n，Wb_n，Wc_n，Wd_n，n＝1，2，...，50的这些系数都分别乘以

计算

其中[·]表示取整运算。将Wa_n，Wb_n，Wc_n，Wd_n，n＝51，52，...，100的这些系数都分别乘以

计算

其中[·]表示取整运算。将Wa_n，Wb_n，Wc_n，Wd_n，n＝101，102，...，150的这些系数都分别乘以

计算

其中[·]表示取整运算。将Wa_n，Wb_n，Wc_n，Wd_n，n＝151，152，...，200的这些系数都分别乘以

计算

其中[·]表示取整运算。将h₀，h₁，h₂，h₃，h₄，h₅，h₆，h₇，h₈，h₉，h₁₀，h₁₁都分别乘以

步骤(d):

将Wa_n，n＝1，2，...，50这些系数按n的顺序放在L₁a₁₀₀₀系数后面，得到L₁a₁，L₁a₂，...，L₁a₁₀₀₀，Wa₁，Wa₂，...，Wa₅₀序列。将Wa_n，n＝51，52，...，100这些系数按n的顺序放在L₂a₁₀₀₀系数后面；将Wa_n，n＝101，102，...，150这些系数按n的顺序放在L₃a₁₀₀₀系数后面；将Wa_n，n＝151，152，...，200这些系数按n的顺序放在L₄a₁₀₀₀系数后面；

将Wb_n，n＝1，2，...，50这些系数按n的顺序放在L₁b₁₀₀₀系数后面，得到L₁b₁，L₁b₂，...，L₁b₁₀₀₀，Wb₁，Wb₂，...，Wb₅₀序列；将Wb_n，n＝51，52，...，100这些系数按n的顺序放在L₂b₁₀₀₀系数后面；将Wb_n，n＝101，102，...，150这些系数按n的顺序放在L₃b₁₀₀₀系数后面；将Wb_n，n＝151，152，...，200这些系数按n的顺序放在L₄b₁₀₀₀系数后面；

同样，使用上述Wa_n和Wb_n分片放置的方法，将Wc_n的系数分片按顺序放到L₁c₁₀₀₀，L₂c₁₀₀₀，L₃c₁₀₀₀，L₄c₁₀₀₀的后面。将Wd_n的系数分片按顺序放到L₁d₁₀₀₀，L₂d₁₀₀₀，L₃d₁₀₀₀，L₄d₁₀₀₀的后面。

将h₀，h₁，h₂，h₃，h₄，h₅，h₆，h₇，h₈，h₉，h₁₀，h₁₁按顺序分别放在L₅a₁₀₀₀，L₅b₁₀₀₀，L₅c₁₀₀₀和L₅d₁₀₀₀的后面，分别得到L₅a₁，L₅a₂，...，L₅a₁₀₀₀，h₀，h₁，...，h₁₀，h₁₁序列；L₅b₁，L₅b₂，...，L₅b₁₀₀₀，h₀，h₁，...，h₁₀，h₁₁序列；L₅c₁，L₅c₂，...，L₅c₁₀₀₀，h₀，h₁，...，h₁₀，h₁₁序列；L₅d₁，L₅d₂，...，L₅d₄₉，L₅d₁₀₀₀，h₀，h₁，...，h₁₀，h₁₁序列；

步骤(e):

对于L₁a₁，L₁a₂，...，L₁a₁₀₀₀，Wa₁，Wa₂，...，Wa₅₀；L₁b₁，L₁b₂，...，L₁b₁₀₀₀，Wb₁，Wb₂，...，Wb₅₀；L₁c₁，L₁c₂，...，L₁c₅₀，Wc₁，Wc₂，...，Wc₁₀₀₀；L₁d₁，L₁d₂，...，L₁d₁₀₀₀，Wd₁，Wd₂，...，Wd₅₀，按公式(5)和公式(6)所示的逆曲线傅里叶变换公式计算顶点坐标：

其中，公式中的a_k按k值由小到大的顺序依次对应了L₁a₁，L₁a₂，...，L₁a₁₀₀₀，Wa₁，Wa₂，...，Wa₅₀。同样，b_k对应了L₁b₁，L₁b₂，...，L₁b₁₀₀₀，Wb₁，Wb₂，...，Wb₅₀；c_k对应了L₁c₁，L₁c₂，...，L₁c₅₀，Wc₁，Wc₂，...，Wc₁₀₀₀；d_k对应了L₁d₁，L₁d₂，...，L₁d₁₀₀₀，Wd₁，Wd₂，...，Wd₅₀。K是系数序列的系数个数。n＝0，1，2，...，N-1，N是L₁曲线的顶点个数。T是TL₁，t₀等于0，t₁是L₁的起点沿着曲线到第一个顶点的曲线长度，t₂是L₁的起点沿着曲线到第二个顶点的曲线长度，后面以此类推，实际上t_n在步骤(b)中已经计算得到。

将公式(5)和公式(6)计算得到的x(n)作为横坐标，y(n)作为纵坐标，得到一个顶点序列(x(n),y(n))，n＝0，1，2，...，N-1。将这个顶点序列依次用线段连接可以得到一个闭合曲线WL₁。把这个闭合曲线WL₁替代L₁保存到等高线地图M中。

对于L₂a₁，L₂a₂，...，L₂a₁₀₀₀，Wa₅₁，Wa₅₂，...，Wa₁₀₀；L₂b₁，L₂b₂，...，L₂b₁₀₀₀，Wb₅₁，Wb₅₂，...，Wb₁₀₀；L₂c₁，L₂c₂，...，L₂c₁₀₀₀，Wc₅₁，Wc₅₂，...，Wc₁₀₀；L₂d₁，L₂d₂，...，L₂d₁₀₀₀，Wd₅₀，Wd₅₁，...，Wd₁₀₀，采用与WL₁同样的方法，利用逆曲线傅里叶变换可以得到一个顶点序列，将这个顶点序列依次连接得到WL₂，将WL₂替代L₂保存到等高线地图M中。

对于L₃a₁，L₃a₂，...，L₃a₁₀₀₀，Wa₁₀₁，Wa₁₀₂，...，Wa₁₅₀；L₃b₁，L₃b₂，...，L₃b₁₀₀₀，Wb₁₀₁，Wb₁₀₂，...，Wb₁₅₀；L₃c₁，L₃c₂，...，L₃c₁₀₀₀，Wc₁₀₁，Wc₁₀₂，...，Wc₁₅₀；L₃d₁，L₃d₂，...，L₃d₁₀₀₀，Wd₁₀₁，Wd₁₀₂，...，Wd₁₅₀，采用与WL₁同样的方法，利用逆曲线傅里叶变换可以得到一个顶点序列，将这个顶点序列依次连接得到WL₃，将WL₃替代L₃保存到等高线地图M中。

对于L₄a₁，L₄a₂，...，L₄a₁₀₀₀，Wa₁₅₁，Wa₁₅₂，...，Wa₂₀₀；L₄b₁，L₄b₂，...，L₄b₅₀，Wb₁₅₁，Wb₁₅₂，...，Wb₂₀₀；L₄c₁，L₄c₂，...，L₄c₁₀₀₀，Wc₁₅₁，Wc₁₅₂，...，Wc₂₀₀；L₄d₁，L₄d₂，...，L₄d₁₀₀₀，Wd₁₅₁，Wd₁₅₂，...，Wd₂₀₀，采用与WL₁同样的方法，利用逆曲线傅里叶变换可以得到一个顶点序列，将这个顶点序列依次连接得到WL₄，将WL₄替代L₄保存到等高线地图M中。

对于，L₅a₁，L₅a₂，...，L₅a₁₀₀₀，h₀，h₁，...，h₁₀，h₁₁；L₅b₁，L₅b₂，...，L₅b₁₀₀₀，h₀，h₁，...，h₁₀，h₁₁；L₅c₁，L₅c₂，...，L₅c₁₀₀₀，h₀，h₁，...，h₁₀，h₁₁；L₅d₁，L₅d₂，...，L₅d₁₀₀₀，h₀，h₁，...，h₁₀，h₁₁，采用与WL₁同样的方法，利用逆曲线傅里叶变换可以得到一个顶点序列，将这个顶点序列依次连接可以得到WL₅，将WL₅替代L₅保存到等高线地图M中。

这样，经过步骤(a)至步骤(e)，将水印信息w＝(w₁，w₂，...，w_m，...，w₁₁)嵌入到了等高线地图M中,得到包含水印信息的等高线地图M′。

(二)、水印信息的提取是从待检测的等高线地图M′中提取水印信息w的内容，输入为包含水印信息的等高线地图M′，输出为等高线地图里包含的水印信息w，具体包括以下步骤：

步骤(1)：

对等高线地图M′的所有闭合曲线分别计算其周长，对所有闭合曲线按周长进行排序，周长越大，排序越靠前。如果有曲线周长相同，那么对周长相等的闭合曲线再按顶点数量进行排序，顶点数量越大，排序越前。

对等高线地图M′的所有闭合曲线排好序后，在排好序的闭合曲线中选择最靠前的五条闭合曲线作为潜入数字水印的对象。这里将排序第一的闭合曲线记为L₁,排序第二的闭合曲线记为L₂，后面以此类推，排序第三，四，五的闭合曲线分别记为L₃，L₄和L₅。将L₁的周长记为TL₁，L₂的周长记为TL₂，L₃的周长记为TL₃，L₄的周长记为TL₄，L₅的周长记为TL₅。

步骤(2)：

步骤(a):

对L₁，从0度水平开始按逆时针顺序扫描闭合曲线所记录包含的顶点，可以获得一个顶点序列。对于该顶点序列，使用公式(1)-公式(4)所示的曲线傅里叶变换公式计算得到一组曲线傅里叶变换系数，n＝1，2，...，1050。这里把公式(1)计算得到的系数记为L₁a_n，把公式(2)计算得到的系数记为L₁b_n，把公式(3)计算得到的系数记为L₁c_n，把公式(4)计算得到的系数记为L₁d_n。

把L₁a₁₀₀₁赋值给Wa₁，把L₁a₁₀₀₂赋值给Wa₂，以此类推，这样得到Wa₁，Wa₂，...，Wa₅₀。把L₁b₁₀₀₁赋值给Wb₁，把L₁b₁₀₀₂赋值给Wb₂，以此类推，这样得到Wb₁，Wb₂，...，Wb₅₀。用同样的方法，利用L₁c_n和L₁d_n分别得到Wc₁，Wc₂，...，Wc₅₀和Wd₁，Wd₂，...，Wd₅₀。

步骤(b):

按如上所述步骤(a)的方法，对于L₂计算得到Wa_n，Wb_n，Wc_n，Wd_n，n＝51，52，...，100；对于L₃计算得到Wa_n，Wb_n，Wc_n，Wd_n，n＝101，102，...，150；对于L₄计算得到Wa_n，Wb_n，Wc_n，Wd_n，n＝151，152，...，200。

步骤(c):

对于L₅,采用步骤(a)L₁的方法得到顶点序列，然后对于该顶点序列，使用公式(1)-公式(4)所示的曲线傅里叶变换公式计算得到一组曲线傅里叶变换系数，n＝1，2，...，1012。最后把L₁a₁₀₀₁赋值给h₀，把L₁a₁₀₀₂赋值给h₁，以此类推，这样得到h₀，h₁，...，h₁₀，h₁₁。

步骤(d):

计算

其中[·]表示取整运算。将Wa_n，Wb_n，Wc_n，Wd_n，n＝1，2，...，50的这些系数都分别除以

计算

其中[·]表示取整运算。将Wa_n，Wb_n，Wc_n，Wd_n，n＝51，52，...，100，这些系数都分别除以

计算

其中[·]表示取整运算。将Wa_n，Wb_n，Wc_n，Wd_n，n＝101，102，...，150，这些系数都分别除以

计算

其中[·]表示取整运算。将Wa_n，Wb_n，Wc_n，Wd_n，n＝151，152，...，200，这些系数都分别除以

计算

其中[·]表示取整运算。将h₀，h₁，...，h₁₀，h₁₁这些系数都分别除以

步骤(3)：

按公式(7)和公式(8)所示的逆曲线傅里叶变换公式计算顶点坐标：

其中，n＝0，1，2，...，11。

计算

n＝1，2，...，11，得到一个实数序列。对实数序列L(n)，n＝1，2，...，11的每个实数分别进行四舍五入计算，然后得到一个Z(n)整数序列。这个整数序列就是提取出来的数字水印信息w。

Claims (4)

1.一种基于曲线傅里叶变换的等高线地图数据盲水印方法，其特征在于该方法包括水印信息嵌入和水印信息提取两部分；

所述的水印信息嵌入是把用户选择的水印信息嵌入到原始等高线地图中得到一个新的等高线地图，输入为一个等高线地图M和一个包含11个0～9之间整数的序列w，输出为包含水印信息的等高线地图M′，具体步骤为：

步骤(1):曲线傅里叶变换计算

设有一顶点序列D₁，D₂，...，D_i，其对应的坐标序列为(x₁,y₁)，(x₂,y₂)，...，(x_i,y_i)，对于顶点序列D₁，D₂，...，D_i，使用曲线傅里叶变换计算得到一组曲线傅里叶变换系数a_n，b_n，c_n和d_n；

步骤(2):调制水印信息

嵌入的水印信息是一个包含11个整数的序列，每个整数取值在1到9之间，记这个整数序列为w＝(w₁，w₂，...，w_m，...，w₁₁)；

画一圆R，其半径为1，在圆R中分别画角度为0度、30度、60度，90度、120度，150度，180度、210度、240度，270度、300度、330度的半径线，这12条半径线与R的交点分别记为A₀、A₁、A₂、A₃、A₄、A₅、A₆、A₇、A₈、A₉、A₁₀及A₁₁；其中，A₀是0度半径线与R的交点，A₁是30度半径线与R的交点，A₃是60度半径线与R的交点，其他以此类推；

把A_k交点在圆R上的坐标分别记为(x_k,y_k),k＝0,1,，...，11那么有：

x_k＝cos(k×π/6)，y_k＝sin(k×π/6)

将w＝(w₁，w₂，...，w_m，...，w₁₁)的11个整数值按从前到后的顺序，分别与A₁、A₂、A₃、A₄、A₅、A₆、A₇、A₈、A₉、A₁₀、A₁₁一一对应，然后把A₁、A₂、A₃、A₄、A₅、A₆、A₇、A₈、A₉、A₁₀、A₁₁的横坐标和纵坐标分别除以w₁，w₂，w₃，w₄，w₅，w₆，w₇，w₈，w₉，w₁₀，w₁₁的整数值；

计算A₀和点A₁之间线段长度h₀；计算A₁和点A₂之间线段长度h₁，计算A₂和A₃之间的线段长度h₂、A₃和A₄之间的线段长度h₃、A₄和A₅之间的线段长度h₄、A₅和A₆之间的线段长度h₅、A₆和A₇之间的线段长度h₆、A₇和A₈之间的线段长度h₇、A₈和A₉之间的线段长度h₈、A₉和A₁₀之间的线段长度h₉、A₁₀和A₁₁之间的线段长度h₁₀及A₁₁和A₀之间的线段长度h₁₁；

对A₀、A₁、A₂、A₃、A₄、A₅、A₆、A₇、A₈、A₉、A₁₀及A₁₁使用曲线傅里叶变换得到一组曲线傅里叶变换系数，得到系数Wa_n、Wb_n、Wc_n和Wd_n；这样，根据水印信息w＝(w₁，w₂，...，w_m，...，w₁₁)，利用曲线傅里叶变换得到一组曲线傅里叶变换系数Wa_n，Wb_n，Wc_n，Wd_n，n＝1，2，...，199，200；

步骤(3)：嵌入水印信息

步骤(a)：

对等高线地图的所有闭合曲线分别计算其周长，对所有闭合曲线按周长进行排序，周长越大，排序越靠前；

对等高线地图的所有闭合曲线排好序后，在排好序的闭合曲线中选择最靠前的五条闭合曲线作为嵌入数字水印的对象；将排序第一的闭合曲线记为L₁,排序第二的闭合曲线记为L₂，排序第三，四，五的闭合曲线分别记为L₃，L₄和L₅；

将L₁的周长记为TL₁，L₂的周长记为TL₂，L₃的周长记为TL₃，L₄的周长记为TL₄，L₅的周长记为TL₅；

步骤(b)：

对L₁，从0度水平开始按逆时针顺序扫描闭合曲线所记录包含的顶点，获得一个顶点序列；对于该顶点序列，使用曲线傅里叶变换得到一组曲线傅里叶变换系数，分别记为L₁a_n、L₁b_n、L₁c_n和L₁d_n，其中n＝1，2，...，999，1000；

在使用曲线傅里叶变换计算L₁顶点序列的曲线傅里叶系数时，同时保存所有的t_j，也就是保存L₁第一个顶点沿着曲线按逆时针方向到其他每一个顶点的曲线距离；

采用与如上L₁同样的计算方式，对于L₂计算得到一组曲线傅里叶变换系数L₂a_n，L₂b_n，L₂c_n，L₂d_n，n＝1，2，...，999，1000和L₂第一个顶点沿着曲线按逆时针方向到其他每一个顶点的曲线距离；

采用与如上L₁同样的方式，对于L₃计算得到一组曲线傅里叶变换系数L₃a_n，L₃b_n，L₃c_n，L₃d_n，n＝1，2，...，999，1000和L₃第一个顶点沿着曲线按逆时针方向到其他每一个顶点的曲线距离；

采用与如上L₁同样的方式，对于L₄计算得到一组曲线傅里叶变换系数L₄a_n，L₄b_n，L₄c_n，L₄d_n，n＝1，2，...，999，1000和L₄第一个顶点沿着曲线按逆时针方向到其他每一个顶点的曲线距离；

采用与如上L₁同样的方式，对于L₅计算得到一组曲线傅里叶变换系数L₅a_n，L₅b_n，L₅c_n，L₅d_n，n＝1，2，...，999，1000和L₅第一个顶点沿着曲线按逆时针方向到其他每一个顶点的曲线距离；

步骤(c):

计算

将Wa_n，Wb_n，Wc_n，Wd_n，n＝1，2，...，50的这些系数都分别乘以

计算

将Wa_n，Wb_n，Wc_n，Wd_n，n＝51，52，...，100的这些系数都分别乘以

计算

将Wa_n，Wb_n，Wc_n，Wd_n，n＝101，102，...，150的这些系数都分别乘以

计算

将Wa_n，Wb_n，Wc_n，Wd_n，n＝151，152，...，200的这些系数都分别乘以

计算

将h₀，h₁，h₂，h₃，h₄，h₅，h₆，h₇，h₈，h₉，h₁₀，h₁₁都分别乘以

其中[·]表示取整运算；

步骤(d):

将Wa_n，n＝1，2，...，50这些系数按n的顺序放在L₁a₁₀₀₀系数后面，得到L₁a₁，L₁a₂，...，L₁a₁₀₀₀，Wa₁，Wa₂，...，Wa₅₀序列；将Wa_n，n＝51，52，...，100这些系数按n的顺序放在L₂a₁₀₀₀系数后面；将Wa_n，n＝101，102，...，150这些系数按n的顺序放在L₃a₁₀₀₀系数后面；将Wa_n，n＝151，152，...，200这些系数按n的顺序放在L₄a₁₀₀₀系数后面；

将Wb_n，n＝1，2，...，50这些系数按n的顺序放在L₁b₁₀₀₀系数后面，得到L₁b₁，L₁b₂，...，L₁b₁₀₀₀，Wb₁，Wb₂，...，Wb₅₀序列；将Wb_n，n＝51，52，...，100这些系数按n的顺序放在L₂b₁₀₀₀系数后面；将Wb_n，n＝101，102，...，150这些系数按n的顺序放在L₃b₁₀₀₀系数后面；将Wb_n，n＝151，152，...，200这些系数按n的顺序放在L₄b₁₀₀₀系数后面；

同样，使用上述Wa_n和Wb_n分片放置的方式，将Wc_n的系数分片按顺序放到L₁c₁₀₀₀，L₂c₁₀₀₀，L₃c₁₀₀₀，L₄c₁₀₀₀的后面；将Wd_n的系数分片按顺序放到L₁d₁₀₀₀，L₂d₁₀₀₀，L₃d₁₀₀₀，L₄d₁₀₀₀的后面；

将h₀，h₁，h₂，h₃，h₄，h₅，h₆，h₇，h₈，h₉，h₁₀，h₁₁按顺序分别放在L₅a₁₀₀₀，L₅b₁₀₀₀，L₅c₁₀₀₀和L₅d₁₀₀₀的后面，分别得到L₅a₁，L₅a₂，...，L₅a₁₀₀₀，h₀，h₁，...，h₁₀，h₁₁序列；L₅b₁，L₅b₂，...，L₅b₁₀₀₀，h₀，h₁，...，h₁₀，h₁₁序列；L₅c₁，L₅c₂，...，L₅c₁₀₀₀，h₀，h₁，...，h₁₀，h₁₁序列；L₅d₁，L₅d₂，...，L₅d₄₉，L₅d₁₀₀₀，h₀，h₁，...，h₁₀，h₁₁序列；

步骤(e):

对于L₁a₁，L₁a₂，...，L₁a₁₀₀₀，Wa₁，Wa₂，...，Wa₅₀；L₁b₁，L₁b₂，...，L₁b₁₀₀₀，Wb₁，Wb₂，...，Wb₅₀；L₁c₁，L₁c₂，...，L₁c₅₀，Wc₁，Wc₂，...，Wc₁₀₀₀；L₁d₁，L₁d₂，...，L₁d₁₀₀₀，Wd₁，Wd₂，...，Wd₅₀，采用逆曲线傅里叶变换得到顶点坐标；

由顶点坐标得到一个顶点序列(x(n),y(n))，n＝0，1，2，...，N-1；将这个顶点序列依次用线段连接得到一个闭合曲线WL₁；把这个闭合曲线WL₁替代L₁保存到等高线地图中；

对于L₂a₁，L₂a₂，...，L₂a₁₀₀₀，Wa₅₁，Wa₅₂，...，Wa₁₀₀；L₂b₁，L₂b₂，...，L₂b₁₀₀₀，Wb₅₁，Wb₅₂，...，Wb₁₀₀；L₂c₁，L₂c₂，...，L₂c₁₀₀₀，Wc₅₁，Wc₅₂，...，Wc₁₀₀；L₂d₁，L₂d₂，...，L₂d₁₀₀₀，Wd₅₀，Wd₅₁，...，Wd₁₀₀，采用与WL₁同样的方式，利用逆曲线傅里叶变换得到一个顶点序列，将这个顶点序列依次连接得到WL₂，将WL₂替代L₂保存到等高线地图中；

对于L₃a₁，L₃a₂，...，L₃a₁₀₀₀，Wa₁₀₁，Wa₁₀₂，...，Wa₁₅₀；L₃b₁，L₃b₂，...，L₃b₁₀₀₀，Wb₁₀₁，Wb₁₀₂，...，Wb₁₅₀；L₃c₁，L₃c₂，...，L₃c₁₀₀₀，Wc₁₀₁，Wc₁₀₂，...，Wc₁₅₀；L₃d₁，L₃d₂，...，L₃d₁₀₀₀，Wd₁₀₁，Wd₁₀₂，...，Wd₁₅₀，采用与WL₁同样的方式，利用逆曲线傅里叶变换得到一个顶点序列，将这个顶点序列依次连接得到WL₃，将WL₃替代L₃保存到等高线地图中；

对于L₄a₁，L₄a₂，...，L₄a₁₀₀₀，Wa₁₅₁，Wa₁₅₂，...，Wa₂₀₀；L₄b₁，L₄b₂，...，L₄b₅₀，Wb₁₅₁，Wb₁₅₂，...，Wb₂₀₀；L₄c₁，L₄c₂，...，L₄c₁₀₀₀，Wc₁₅₁，Wc₁₅₂，...，Wc₂₀₀；L₄d₁，L₄d₂，...，L₄d₁₀₀₀，Wd₁₅₁，Wd₁₅₂，...，Wd₂₀₀，采用与WL₁同样的方式，利用逆曲线傅里叶变换得到一个顶点序列，将这个顶点序列依次连接得到WL₄，将WL₄替代L₄保存到等高线地图中；

对于，L₅a₁，L₅a₂，...，L₅a₁₀₀₀，h₀，h₁，...，h₁₀，h₁₁；L₅b₁，L₅b₂，...，L₅b₁₀₀₀，h₀，h₁，...，h₁₀，h₁₁；L₅c₁，L₅c₂，...，L₅c₁₀₀₀，h₀，h₁，...，h₁₀，h₁₁；L₅d₁，L₅d₂，...，L₅d₁₀₀₀，h₀，h₁，...，h₁₀，h₁₁，采用与WL₁同样的方式，利用逆曲线傅里叶变换得到一个顶点序列，将这个顶点序列依次连接可以得到WL₅，将WL₅替代L₅保存到等高线地图中；

这样，经过步骤(a)至步骤(e)，将水印信息w＝(w₁，w₂，...，w_m，...，w₁₁)嵌入到了等高线地图中；

所述的水印信息的提取是从待检测的等高线地图M′中提取水印信息w的内容，输入为包含水印信息的等高线地图M′，输出为等高线地图里包含的水印信息w，具体包括以下步骤：

步骤(1)：

对等高线地图的所有闭合曲线分别计算其周长，对所有闭合曲线按周长进行排序，周长越大，排序越靠前；

对等高线地图的所有闭合曲线排好序后，在排好序的闭合曲线中选择最靠前的五条闭合曲线作为潜入数字水印的对象；将排序第一的闭合曲线记为L₁,排序第二的闭合曲线记为L₂，排序第三，四，五的闭合曲线分别记为L₃，L₄和L₅；将L₁的周长记为TL₁，L₂的周长记为TL₂，L₃的周长记为TL₃，L₄的周长记为TL₄，L₅的周长记为TL₅；

步骤(2)：

步骤(a):

对L₁，从0度水平开始按逆时针顺序扫描闭合曲线所记录包含的顶点，获得一个顶点序列；对于该顶点序列，使用曲线傅里叶变换得到一组曲线傅里叶变换系数，记为L₁a_n、L₁b_n、L₁c_n和L₁d_n，n＝1，2，...，1050；

把L₁a₁₀₀₁赋值给Wa₁，把L₁a₁₀₀₂赋值给Wa₂，以此类推，这样得到Wa₁，Wa₂，...，Wa₅₀；把L₁b₁₀₀₁赋值给Wb₁，把L₁b₁₀₀₂赋值给Wb₂，以此类推，这样得到Wb₁，Wb₂，...，Wb₅₀；用同样的方式，利用L₁c_n和L₁d_n分别得到Wc₁，Wc₂，...，Wc₅₀和Wd₁，Wd₂，...，Wd₅₀；

步骤(b):

按如上所述步骤(a)的方式，对于L₂计算得到Wa_n，Wb_n，Wc_n，Wd_n，n＝51，52，...，100；对于L₃计算得到Wa_n，Wb_n，Wc_n，Wd_n，n＝101，102，...，150；对于L₄计算得到Wa_n，Wb_n，Wc_n，Wd_n，n＝151，152，...，200；

步骤(c):

对于L₅,采用步骤(a)L₁的方式得到顶点序列，然后对于该顶点序列，使用曲线傅里叶变换公式计算得到一组曲线傅里叶变换系数，n＝1，2，...，1012；最后把L₁a₁₀₀₁赋值给h₀，把L₁a₁₀₀₂赋值给h₁，以此类推，这样得到h₀，h₁，...，h₁₀，h₁₁；

步骤(d):

计算

将Wa_n，Wb_n，Wc_n，Wd_n，n＝1，2，...，50的这些系数都分别除以

计算

将Wa_n，Wb_n，Wc_n，Wd_n，n＝51，52，...，100，这些系数都分别除以

计算

将Wa_n，Wb_n，Wc_n，Wd_n，n＝101，102，...，150，这些系数都分别除以

计算

将Wa_n，Wb_n，Wc_n，Wd_n，n＝151，152，...，200，这些系数都分别除以

计算

将h₀，h₁，...，h₁₀，h₁₁这些系数都分别除以

步骤(3)：

按如下公式计算顶点坐标：

其中，n＝0，1，2，...，11；

计算

n＝1，2，...，11，得到一个实数序列；对实数序列L(n)，n＝1，2，...，11的每个实数分别进行四舍五入计算，然后得到一个Z(n)整数序列；这个整数序列就是提取出来的数字水印信息。

2.根据权利要求1所述的一种基于曲线傅里叶变换的等高线地图数据盲水印方法，其特征在于：曲线傅里叶变换计算采用如下公式(1)到公式(4)

其中，Δx_j＝x_j-x_j-1，Δy_j＝y_j-y_j-1，

T＝t_i，

n＝1,2,3...。

3.根据权利要求1所述的一种基于曲线傅里叶变换的等高线地图数据盲水印方法，其特征在于：在嵌入水印信息的步骤(a)中，如果有曲线周长相同，那么对周长相等的闭合曲线再按顶点数量进行排序，顶点数量越大，排序越前。

4.根据权利要求1所述的一种基于曲线傅里叶变换的等高线地图数据盲水印方法，其特征在于：在嵌入水印信息的步骤(e)中逆曲线傅里叶变换采用如下公式：

其中，公式中的a_k按k值由小到大的顺序依次对应L₁a₁，L₁a₂，...，L₁a₁₀₀₀，Wa₁，Wa₂，...，Wa₅₀；同样，b_k对应L₁b₁，L₁b₂，...，L₁b₁₀₀₀，Wb₁，Wb₂，...，Wb₅₀；c_k对应L₁c₁，L₁c₂，...，L₁c₅₀，Wc₁，Wc₂，...，Wc₁₀₀₀；d_k对应L₁d₁，L₁d₂，...，L₁d₁₀₀₀，Wd₁，Wd₂，...，Wd₅₀；K是系数序列的系数个数；n＝0，1，2，...，N-1，N是L₁曲线的顶点个数；T是TL₁，t₀等于0，t₁是L₁的起点沿着曲线到第一个顶点的曲线长度，t₂是L₁的起点沿着曲线到第二个顶点的曲线长度，以此类推。