太牛了哥，我想知道您的这些AIS数据是从哪里下载的？是需要收费的吗？

一楼附件里面有测试数据的

轨迹质量提升策略

由于传输媒介和设备等因素的影响，AIS信号可能会发生中断和丢失。这导致在时空维度上，相邻位置点之间存在过大的采样间隔，使得它们彼此相距甚远。此外，原始轨迹数据集缺乏船舶的航次信息，如果不对具有相同MMSI编号的船舶轨迹进行分段，可能会产生错误的轨迹聚类结果。

1. 轨迹分段方法

   轨迹分段的目的仅是为轨迹去噪提供必要的前提。为了对原始轨迹进行分段，将时间阈值设置为2小时，距离阈值设置为2海里。轨迹分段结果可能会产生一些轨迹点数量较少或轨迹长度过短的轨迹段（少于5个轨迹点的轨迹段被视为噪声），这将给后续的去噪、压缩带来不必要的干扰，因此，这些轨迹段应该被视为噪声并且进行滤除。以原始数据voyage_id=106998为例，使用以下代码进行分段（图2.4），其中voyagePartition函数可从附件下载（见1楼），见Segmentation.py。

   

   （图2.4以voyage_id=106998为例的分段代码）

   图2.5为使用两种颜色表示分段后的轨迹线。

   

   

   （图2.5 分段前后对比图） 

2.  轨迹去噪方法
   

   为了滤除海量AIS数据中存在的单个或多个连续的轨迹噪声点，基于Anomaly Detection and Restoration for AIS Raw Data提到的轨迹修复算法，利用船舶自身轨迹信息进行异常点判断与修复的方法，该方法综合考虑经纬度、速度、加速度和航向数据进行异常点判定，相比基于单一位置数据的判定方法，该方法能够有效减少漏判，并且无需依赖历史数据，适用范围更广。我们主要复现了论文中的“异常加速点”和“异常漂移点”的部分，有兴趣的读者可以仔细研读该文章，它还涉及异常转向点和轨迹插值，然后复现这两部分。它的核心封装代码如下，完整代码见附件Clean.py，去噪代码如图2.6所示。

   # 1. 轨迹点之间的haversine距离
   def haversine(self, lon1, lat1, lon2, lat2):
           """
               Calculate the great circle distance between two points
               on the earth (specified in decimal degrees)
           """
           # convert decimal degrees to radians
           lon1, lat1, lon2, lat2 = map(radians, [lon1, lat1, lon2, lat2])
           # haversine formula
           dlon = lon2 - lon1
           dlat = lat2 - lat1
           a = sin(dlat / 2) ** 2 + cos(lat1) * cos(lat2) * sin(dlon / 2) ** 2
           c = 2 * asin(sqrt(a))
           km = 6371 * c
           return km * 1000
   
   # 2. 初始条件
   def __init__(self, df: pd.DataFrame, vd: float = 25):
           self.traj = df
           self.traj_length = len(self.traj)
           # vd: ship design speed, unit: m/s
           self.vd = vd * (1852 / 3600)
           # The ship length is used to calculate the maximum and minimum acceleration
           self.length = self.traj["length"].unique()[0]
           self.max_length = 2 * 10 * self.length
           self.min_length = 2 * 8 * self.length
           # max_acc、min_acc: maximum and minimum acceleration of the ship, unit: m/s^2
           self.max_acc = self.vd ** 2 / self.max_length
           self.min_acc = -self.vd ** 2 / self.min_length
   
   # 3. 异常轨迹点处理
   def _coordinatesAbnormal(self, traj_df):
           traj_df["timestamp"] = traj_df["updatetime"].apply(lambda x: mktime(x.timetuple()))
           toDeleteIndex = set()
           for i, rows in traj_df.iterrows():
               if i == len(traj_df) - 1:
                   break
               nextRow = traj_df.iloc[i+1]
               vi, vj = rows["speed"] * (1852 / 3600), nextRow["speed"] * (1852 / 3600)
               ti, tj = rows.timestamp, nextRow.timestamp
               try:
                   tm = (vj - vi - self.min_acc * tj + self.max_acc * ti) / (self.max_acc - self.min_acc)
                   smax = vi * (tm - ti) + vj * (tj - tm) + \
                          0.5 * self.max_acc * pow(tm - ti, 2) - \
                          0.5 * self.min_acc * pow(tj - tm, 2)
                   # distance, unit: m
                   distance = self.haversine(rows.lon, rows.lat, nextRow.lon, nextRow.lat)
                   if distance > smax:
                       toDeleteIndex.add(i + 1)
               except ZeroDivisionError as e:
                   print("acc_max is equal to acc_min")
           coords_df = traj_df.drop(toDeleteIndex).reset_index(drop=True)
           return coords_df.drop("timestamp", axis=1)
   
   # 4. 异常航速点处理
   def _speedAbnormal(self, traj_df):
           toDeleteIndex = set()
           for i, rows in traj_df.iterrows():
               if i == len(traj_df) - 2:
                   break
               nextRow = traj_df.iloc[i + 1]
               deltaSpeed = (nextRow["speed"] - rows["speed"]) * (1852 / 3600)
               dt = (nextRow["updatetime"] - rows["updatetime"]).total_seconds()
               acc = deltaSpeed / dt
               current_speed = rows["speed"]
               next_speed = nextRow["speed"]
               current_before_speed = traj_df.iloc[i - 1]["speed"]
               next_next_speed = traj_df.iloc[i + 2]["speed"]
   
               if (self.max_acc - acc) * (acc - self.min_acc) < 0:
                   if abs(current_speed - current_before_speed) > abs(next_speed - next_next_speed):
                       toDeleteIndex.add(i)
                   else:
                       toDeleteIndex.add(i + 1)
           return traj_df.drop(toDeleteIndex).reset_index(drop=True)

   

   （图2.6 轨迹去噪代码）

   去噪后，轨迹点数量变为581,879。去噪前后对比图2.7如下所示：

   

   

   （图2.7 轨迹去噪前后对比）

3. 轨迹压缩方法
   

   轨迹压缩的主要目标是找到轨迹曲线上的关键特征点，一般用于研究水域转向点的提取，结合点轨迹聚类算法，也可识别出研究水域内的转向区域。压缩算法的另一个应用是使用这些轨迹特征点近似地表达轨迹的原始形状，减轻后续轨迹聚类中轨迹相似度的计算开销，大大增加轨迹聚类效率。我们将使用Douglas–Peucker（道格拉斯-普克算法）来对船舶轨迹进行压缩，图2.8是道格拉斯-普克算法的简要描述，它的核心处理过程如下：

   （1）连接轨迹曲线首尾两点A，B
   （2）依次计算曲线上所有点到直线AB的距离
   （3）计算得到最大距离dmax。如果dmax小于给定压缩阈值，则去除曲线上除A、B的所有点；反之，以最大距离对应的点为界，将轨迹曲线分为两段
   （4）对轨迹曲线分段重复1）-3）步骤，直到所有dmax均小于给定压缩阈值，即完成轨迹压缩过程

   

   （图2.8 轨迹压缩过程）

   这里我们将给出实现轨迹压缩的完整代码供读者学习。

   # 1. 求点到直线的距离，见上述步骤（1）-（3）,用到第三方库numpy和Python的自带数学计算库math（import numpy和from math import *）
   def pointLineDistance(self, point, linePoint1, linePoint2):
           """
               点到直线的距离
           """
           x0, y0 = point[0], point[1]
           x1, y1 = linePoint1[0], linePoint1[1]
           x2, y2 = linePoint2[0], linePoint2[1]
           numerator = abs((x2 - x1) * (y1 - y0) - (x1 - x0) * (y2 - y1))
           denominator = sqrt(pow((x2 - x1), 2) + pow((y2 - y1), 2))
           distance = numerator / denominator
           return distance
   
   # 2. 递归完成轨迹压缩
   def douglasPeucker(self, pointList, epsilon):
           """
               Douglas-Peucker
           """
           dmax = 0
           index = 0
           end = len(pointList)
           for i in range(1, end-1):
               d = self.pointLineDistance(pointList[i][:2], pointList[0][:2], pointList[-1][:2])
               if d > dmax:
                   dmax = d
                   index = i
   
           if dmax > epsilon:
               recursive1 = self.douglasPeucker(pointList[: index + 1], epsilon)
               recursive2 = self.douglasPeucker(pointList[index: ], epsilon)
               return np.vstack((recursive1[:-1], recursive2))
           else:
               return np.stack((pointList[0], pointList[-1]))

   完整代码见附件douglasPeucker.py。实验数据可从原始数据中筛选出voyage_id=106998的轨迹，自己调整epsilon参数，实现不同压缩参数的轨迹压缩，结合QGIS可视化工具，观察效果。

大佬我有看到您在1l里的数据，但我想知道您这个数据是从哪里获得的？我在这个网站里没有找的？是需要购买还是有部分数据是免费开源的

其它数据是需要购买的。如需要，请联系客服