Method for vessel traffic pattern recognition via data quality control and data compression

Abstract

The present invention proposes a vessel traffic pattern identification method via data quality control and data compression. Firstly, assort a collection of Automatic Identification system (AIS) data points according to Mobile Service Identify (MMSI) code and sort each collection result by time ascending order, and then delete duplicated vessel AIS data points considering time stamp, latitude, longitude and vessel speed over ground, then segment vessel trajectories. Secondly obtain high-quality AIS data with an AIS data anomaly detection and repair and compress each vessel trajectory with the Douglas-Peucker algorithm. Thirdly, cluster vessel trajectories with the Quick Bundles algorithm, and identify maritime traffic pattern. The invention can efficiently identify vessel traffic patterns, and help maritime traffic management departments to accurately identify a traffic situation.

Claims

1. A method for vessel traffic pattern recognition via data quality control and data compression, comprising the following steps: (1) assorting a collection of AIS data points according to MMSI and sorting each collection result by time ascending order, deleting duplicative AIS data points and segmenting vessel trajectories: allocating each AIS data point in a collection to a vessel trajectory trajectory.sub.z so that each point therein having a same MMSI, and sorting each vessel trajectory trajectory.sub.z by time ascending order, thus obtaining a set of vessel trajectories trajectory = {trajectory.sub.z}, z = 1,2,3, ...,w, wherein trajectory.sub.z denoting a zth vessel trajectory which z = 1,2,3, ...,w, each AIS data point of a vessel trajectory trajectory.sub.z represented by e = {MMSI, Time, Ion, lat, sog} , MMSI denoting a Maritime Mobile Service Identify of vessel, Time denoting a time stamp, Ion denoting a longitude, lat denoting a latitude, and sog denoting a vessel speed over ground for said each vessel trajectory trajectory.sub.z; deleting duplicative AIS data points and segmenting vessel trajectory for each vessel trajectory trajectory.sub.z as follows: for AIS data points therein having a same time stamp, a same longitude, a same latitude, and a same vessel speed over ground, retaining only one thereof, while deleting the others thereof; thereafter segmenting the vessel trajectory trajectory.sub.z: starting from index 1 in trajectory.sub.z to obtain a first AIS data point efirst(j — 1) and a last AIS data point elast(j) such that AIS data points therebetween satisfying Expression set (1), continuing till end of index of trajectory.sub.z while deleting all the AIS data points between efirst(j — 1) and elast(j), segmenting the vessel trajectory trajectory.sub.z at the last AIS data point elast(j); obtaining a new set of vessel trajectories tra = {tra.sub.i}, i = 1,2,3, ... n, wherein tra.sub.i denoting an ith vessel trajectory with each AIS data point of the vessel trajectory tra.sub.i represented by e = {MMSI, Time, Ion, lat, sog}; $\left\{\begin{array}{c}{\text{sog}}_{\text{j}}<1\\ {\text{time}}_{\text{elast}\left(\text{j}\right)}-{\text{time}}_{\text{efirst}\left(\text{j}-\text{1}\right)}>{\text{Time}}_{\max }\end{array}\right)$ wherein sog.sub.j denoting a speed over ground at a jth AIS data point in the vessel trajectory trajectory.sub.z, time.sub.efirst(j-1) denoting a timestamp of an AIS data point efirst(j — 1) in the vessel trajectory trajectory.sub.z, time.sub.elast(.sub.j) denoting a timestamp of an AIS data point elast(j) in the vessel trajectory trajectory.sub.z, and Time.sub.max denoting a pre-set time threshold; (2) identifying adrift AIS data points and missing vessel trajectory segments for each vessel trajectory tra.sub.i, repairing the missing vessel trajectory segments with cubic spline interpolation algorithm after deleting the adrift AIS data points for said each vessel trajectory tra.sub.i as follows: (2.1) deleting an adrift AIS data point e.sub.j which satisfying Expression set (2): $\left\{\begin{array}{c}\Delta {\text{t}}_{\text{j}}={\text{Time}}_{\text{j}}-{\text{Time}}_{\text{j-1}}\\ \Delta {\text{d}}_{\text{j}}={\text{speed}}_{\max }\ast \Delta {\text{t}}_{\text{j}}\\ \Delta {\text{lon}}_{\text{j}}={\text{lon}}_{\text{j}}-{\text{lon}}_{\text{j}-\text{1}}\ge \frac{\Delta {\text{d}}_{\text{j}}}{\cos {30}^{\circ }\ast \frac{\text{π}}{180}\ast 111000}\\ \Delta {\text{lat}}_{\text{j}}={\text{lat}}_{\text{j}}-{\text{lat}}_{\text{j}-1}\ge \frac{\Delta {\text{d}}_{\text{j}}}{111000}\\ \Delta {\text{t}}_{\text{j+1}}={\text{Time}}_{\text{j+1}}-{\text{Time}}_{\text{j}}\\ \Delta {\text{d}}_{\text{j+1}}={\text{Speed}}_{\max }\ast \Delta {\text{t}}_{\text{j+1}}\\ \Delta {\text{lon}}_{\text{j+1}}={\text{lon}}_{\text{j+1}}-{\text{lon}}_{\text{j}}\ge \frac{\Delta {\text{d}}_{\text{j+1}}}{\cos {30}^{\circ }\ast \frac{\text{π}}{180}\ast 111000}\\ \Delta {\text{lat}}_{\text{j+1}}={\text{lat}}_{\text{j+1}}-{\text{lat}}_{\text{j}}\ge \frac{\Delta {\text{d}}_{\text{j+1}}}{111000}\end{array}\right)$ wherein Δt.sub.j denoting a time interval from adjacent AIS data points e.sub.j-1 to e.sub.j in a vessel trajectory, Time.sub.j-1 denoting a time stamp of an AIS data point e.sub.j-1, Time.sub.j denoting a time stamp of an AIS data point e.sub.j, Δt.sub.j+1 denoting a time interval from adjacent AIS data points e.sub.j+1 to e.sub.j in a vessel trajectory, Time.sub.j+1 denoting a time stamp of an AIS data point e.sub.j+1; (2.2) identifying missing vessel trajectory segments with Expression set (3) wherein a time interval Δt between adjacent AIS data points being greater than 3 min and less than 5 min; $\left\{\begin{array}{c}\Delta \text{t}\text{=}{\text{Time}}_{\text{j+1}}-{\text{Time}}_{\text{j}}\\ 3\min <\Delta \text{t}\text{<}\text{5}\text{min}\end{array}\right)$ (2.3) repairing the missing vessel trajectory segments by cubic spline interpolation algorithm in Eq. (4) subsequent to deletion of the adrift AIS data points in step (2.1) to obtain high-quality AIS data, for each missing vessel trajectory segment as follows: dividing a time series [A, B] of missing vessel trajectory segment into u intervals according to a time interval of 30 seconds, namely [[x.sub.1, x.sub.2], [x.sub.2, x.sub.3], ..., [x.sub.u, x.sub.u+1]] , each sub-time series [x.sub.1, x.sub.2], [x.sub.2, x.sub.3], ..., [x.sub.u-1, x.sub.u] with 30 seconds time interval, a time interval of a sub-time series [x.sub.u, x.sub.u+1] being less than or equal to 30 seconds, A ≤ x.sub.1 < x.sub.2 < ... < x.sub.u < x.sub.u+1 ≤ B; x.sub.1,x.sub.2,x.sub.3, ...,x.sub.u+1 corresponding to function values of y.sub.1,y.sub.2,y.sub.3, .sub....,y.sub.u+1 with y.sub.U = S(x.sub.U), (U = 1,2, ...,u), each sub-time series [x.sub.U, x.sub.U+1] satisfying Eq. (4); interpolating a longitude Ion and a latitude lat and a vessel speed over ground sog of each time point x.sub.U in the missing vessel trajectory segment, y denoting a longitude Ion when interpolating a longitude of a time point, y denoting a latitude lat when interpolating a latitude of a time point, y denoting a vessel speed over ground sog when interpolating a vessel speed over ground of a time point, obtaining a new vessel track.sub.i after a vessel trajectory repair; ${\text{S}}_{\text{U}}\left(\text{x}\right)\text{}\text{=}{\text{a}}_{\text{U}}{\text{x}}^{\text{3}}+{\text{b}}_{\text{U}}{\text{x}}^{\text{2}}+{\text{c}}_{\text{U}}\text{x}\text{+}{\text{d}}_{\text{U}}$ wherein a.sub.U, b.sub.U, c.sub.U, d.sub.U denoting pending coefficients which being derived from the missing vessel trajectory segment; obtaining a new set of vessel trajectories track = {track.sub.i}, i = 1,2,3, ... n after processing each vessel trajectories tra.sub.i in step (2), wherein track.sub.i denoting a ith vessel trajectory in track which i = 1,2,3, ... n, each AIS data point of a vessel trajectory track.sub.i represented by e = {MMSI, Time, Ion, lat, sog}; (3) compressing each vessel trajectory track.sub.i with a Douglas-Peucker algorithm by means of a self-invoking computer program as step (3.3) as follows: (3.1) forming a set of vessel trajectory points p = {p.sub.j(lon.sub.j,lat.sub.j)},j = 1,2,3, ...,v from the vessel trajectory track.sub.i, wherein p.sub.j denoting a jth vessel trajectory point for j = 1,2,3, ...,v, lon.sub.j denoting a jth longitude value in vessel trajectory point p.sub.j, lat.sub.j denoting a jth latitude value in vessel trajectory point p.sub.j; converting each vessel trajectory point p.sub.j from longitude and latitude coordinates to a Mercator coordinates vessel trajectory point m.sub.j with Equation set (5), thus obtaining M = {m.sub.j(mlon.sub.j,mlat.sub.j)},j = 1,2,3, ..., v, wherein M denoting a set of vessel trajectory points in the Mercator coordinate system and M = {m.sub.1 (mlon.sub.1, mlat.sub.1), m.sub.2 (mlon.sub.2, mlat.sub.2), m3 (mlon.sub.3, mlat.sub.3), ..., m.sub.v(mlon.sub.v, mlat.sub.v)} , m.sub.j denoting a jth vessel trajectory point in the Mercator coordinate system which j = 1,2,3, ...,v, mlon.sub.j denoting a jth longitude value in vessel trajectory point m.sub.j in Mercator coordinate system, mlat.sub.j denoting a jth latitude value in vessel trajectory point m.sub.j in the Mercator coordinate system; $\left\{\begin{array}{c}\text{radius}=\frac{\text{lr}\ast \cos \beta }{\sqrt{1-{\text{E}}^2\ast {\sin }^2\beta }}\\ {\text{q}}_{\text{j}}=\ln \left(\tan \left(\frac{\text{π}}{4}+\frac{{\text{lat}}_{\text{j}}}{2}\right){\left(\frac{1-\text{E}\ast \sin {\text{lat}}_{\text{j}}}{1+\text{E}\ast \sin {\text{lat}}_{\text{j}}}\right)}^2\right)\\ {\text{Mlon}}_{\text{j}}=\text{radius}\ast {\text{lon}}_{\text{j}}\\ {\text{Mlat}}_{\text{j}}=\text{radius}\ast {\text{q}}_{\text{j}}\end{array}\right)$ wherein radius denoting a radius of the standard latitude-parallel circle, lr denoting a long radius of Earth’s ellipsoid, β a standard latitude in the Mercator projection, E denoting a first eccentricity of Earth’s ellipsoid, q.sub.j denoting an equivalent latitude of a jth vessel trajectory point; (3.2) initiating in respective of the set of vessel trajectory points M = {m.sub.1(mlon.sub.1, mlat.sub.1), m.sub.2(mlon.sub.2, mlat.sub.2), m.sub.3(mlon.sub.3, mlat.sub.3), ..., m.sub.v(mlon.sub.v, mlat.sub.v)} as follows: denoting r as a set of key vessel trajectory points, putting a starting vessel trajectory point m.sub.1(mlon.sub.1, mlat.sub.1) and an end vessel trajectory point m.sub.v(mlon.sub.v, mlat.sub.v) in the set of vessel trajectory points M as key vessel trajectory points to the set of key vessel trajectory points r in order, obtaining r = {m.sub.1(mlon.sub.1, mlat.sub.1), m.sub.v(mlon.sub.v, mlat.sub.v)}; connecting the starting vessel trajectory point m.sub.1(mlon.sub.1, mlat.sub.1) and the end vessel trajectory point m.sub.v(mlon.sub.v, mlat.sub.v) in the set of vessel trajectory points M as a straight line l.sub.1v, calculating distances dist = {dist.sub.2, dist.sub.3, ..., dist.sub.v-1} from all vessel trajectory points between m.sub.1(mlon.sub.1, mlat.sub.1) and m.sub.v(mlon.sub.v, mlat.sub.v) to the straight line l.sub.1v with Eq. (6), determining a vessel trajectory point m.sub.g(mlon.sub.g, mlat.sub.g) such that dist.sub.g = max {dist.sub.2, dist.sub.3, ..., dist.sub.v-1}; $\text{dist}=\frac{\left|\text{se}\ast \text{ta}\right|}{\left|\text{se}\right|}$ wherein dist denoting a vertical distance from a vessel trajectory point to a straight line in the Mercator coordinate system, se denoting a vector from a start of the straight line to an end of the straight line, ta denoting a vector from the start of the straight line to a target point; concluding step(3.2) on condition dist.sub.g being less than a set compression threshold θ; otherwise, putting the vessel trajectory point m.sub.g(mlon.sub.g, mlat.sub.g) as a key vessel trajectory point to r in order, obtaining r = {m.sub.1 (mlon.sub.1, mlat.sub.1), m.sub.g(mlon.sub.g, mlat.sub.g), m.sub.v(mlon.sub.v, mlat.sub.v)}, dividing the set of vessel trajectory points M = {m.sub.1(mlon.sub.1, mlat.sub.1), m.sub.2(mlon.sub.2, mlat.sub.2), m.sub.3(mlon.sub.3, mlat.sub.3), ..., m.sub.v(mlon.sub.v, mlat.sub.v)} into two sub vessel trajectory point sets M.sub.gsub.sub.h, h = 1,2 from m.sub.1(mlon.sub.1, mlat.sub.1) to m.sub.g(mlon.sub.g, mlat.sub.g) and from m.sub.g(mlon.sub.g, mlat.sub.g) to m.sub.v(mlon.sub.v, mlat.sub.v) , M.sub.gsub.sub.1 = {m.sub.1(mlon.sub.1, mlat.sub.1), ..., m.sub.g(mlon.sub.g, mlat.sub.g)} and M.sub.gsub.sub.2 = {m.sub.g(mlon.sub.g,mlat.sub.g), ...,m.sub.v(mlon.sub.v,mlat.sub.v)}, wherein M.sub.gsub.sub.1 denoting a first set of sub vessel trajectory points, M.sub.gsub.sub.2 denoting a 2nd set of sub vessel trajectory points; calculating a number of vessel trajectory points M.sub.gsub.sub.1number.sub.1 in M.sub.gsub.sub.1 and a number of vessel trajectory points M.sub.gsub.sub.1number.sub.2 in M.sub.gsub.sub.2, processing M.sub.gsub.sub.1 by step (3.3) if the number of vessel trajectory points M.sub.gsub.sub.1number.sub.1 being greater than a set number threshold .Math.; processing M.sub.gsub.sub.2 by step (3.3) if the number of vessel trajectory points M.sub.gsub.sub.1number.sub.2 being greater than the set number threshold .Math.; (3.3) Mtrack = {m.sub.start(mlon.sub.start, mlat.sub.start), ..., m.sub.end(mlon.sub.end, mlat.sub.end)} denoting a sub vessel trajectory point set, m.sub.start(mlon.sub.start, mlat.sub.start) denoting a first vessel trajectory point which start = 1,2,3, ...,v — 1, m.sub.end(mlon.sub.end, mlat.sub.end) denoting a last vessel trajectory point which end = 2,3, ...,v, a subscript start being less than subscript point end; connecting the first point m.sub.start(mlon.sub.start, mlat.sub.start) and the last point m.sub.end(mlon.sub.end, mlat.sub.end) as a straight line l.sub.startend, calculating distances dist = {dist.sub.start+1,dist.sub.start+2,..., dist.sub.end-1,} from all vessel trajectory points between m.sub.start(mlon.sub.start, mlat.sub.start) and m.sub.end(mlon.sub.end, mlat.sub.end) to the straight line l.sub.startend with Eq. (6), determining a vessel trajectory point m.sub.d(mlon.sub.d, mlat.sub.d) such that dist.sub.d = max{dist.sub.start+1 dist.sub.start+2, ..., dist.sub.end-1}, concluding step (3.3) on condition dist.sub.d being less than the compression threshold θ; otherwise, putting the vessel trajectory point m.sub.d(mlon.sub.d, mlat.sub.d) as a key vessel trajectory point to r, dividing the sub vessel trajectory point set Mtrack into two sub vessel trajectory point sets M.sub.dsub.sub.h, h = 1,2 from m.sub.start(mlon.sub.start, mlat.sub.start) to m.sub.d(mlon.sub.d, mlat.sub.d) and m.sub.d(mlon.sub.d, mlat.sub.d) to m.sub.end(mlon.sub.end, mlat.sub.end), M.sub.dsub.sub.1 = {m.sub.start(mlon.sub.start, mlat.sub.start), ..., m.sub.d(mlon.sub.d, mlat.sub.d)} and M.sub.dsub.sub.2 = {m.sub.d(mlon.sub.d, mlat.sub.d), ..., m.sub.end(mlon.sub.end, mlat.sub.end)}, wherein M.sub.dsub.sub.1 denoting a first set of sub vessel trajectory points after splitting the sub vessel trajectory point set Mtrack with the vessel trajectory point m.sub.d(mlon.sub.d, mlat.sub.d) as a split point, M.sub.dsub.sub.2 denoting a 2nd set of sub vessel trajectory points after splitting the sub vessel trajectory point set Mtrack with the vessel trajectory point m.sub.d(mlon.sub.d, mlat.sub.d) as a split point; calculating a number of vessel trajectory points M.sub.dsub.sub.1number.sub.1 in M.sub.dsub.sub.1 and a number of vessel trajectory points M.sub.dsub.sub.1number.sub.2 in M.sub.dsub.sub.2 , processing M.sub.dsub.sub.1 by step (3.3) if the number of vessel trajectory points M.sub.dsub.sub.1number.sub.1 being greater than a set number threshold .Math., processing M.sub.dsub.sub.2 by step (3.3) if the number of vessel trajectory points M.sub.dsub.sub.1number.sub.2 being greater than the set number threshold .Math. until the subscript start greater being than or equal to end; obtaining a new set of vessel trajectories R = {r.sub.i}, i = 1,2,3, ...n after processing each vessel trajectory track.sub.i in step (3), wherein r.sub.i denoting a vessel trajectory of ith vessel which i = 1,2,3, ... n , each vessel trajectory points of vessel trajectory r.sub.i represented by m = {mlon, mlat}; (4) reconstructing each vessel trajectory r.sub.i with cubic spline interpolation algorithm, and clustering vessel trajectories into various clusters by Quick Bundles algorithm to form a vessel traffic pattern as follows: (4.1) reconstructing each vessel trajectory r.sub.i with cubic spline interpolation algorithm, for each vessel trajectory r.sub.i in R, searching a vessel trajectory r.sub.j with most vessel trajectory points, calculating number differences between remaining vessel trajectories and the vessel trajectory r.sub.j trajectory points respectively, and interpolating at an end of each remaining vessel trajectory with cubic spline interpolation algorithm so that each vessel trajectory therein having a same number of trajectory points, obtaining a new set of vessel trajectories T = {T.sub.i{t.sub.j(mlon.sub.j, mlat.sub.j)|j = 1,2,3, ..., k}},i = 1,2,3, ...n, wherein T.sub.i denoting an i th vessel trajectory which i = 1,2,3, ... n, each vessel trajectory T.sub.i being a K × 2 matrix; t.sub.j denoting an j th vessel trajectory point of time order serial number j = 1,2,3, ..., k, each vessel trajectory point t.sub.j of a vessel trajectory T.sub.i represented by t = {mlon, mlat}; each vessel trajectory T.sub.i = (t.sub.1, t.sub.2, •••, t.sub.K) has two ordered polylines, namely a isotropic trajectory T.sub.i = (t.sub.1, t.sub.2, ••• t.sub.K) and a reverse trajectory flip version T.sub.Fi = (t.sub.K, t.sub.K-1, ••• t.sub.1); (4.2) clustering vessel trajectory T.sub.i into various clusters by Quick Bundles algorithm to form a vessel traffic pattern: constructing a cluster class set of vessel trajectories C = {c.sub.q(I, h, s)|q = 1,2, ..., W}, wherein c.sub.q denoting a cluster set of vessel trajectories in cluster q which q = 1,2, ..., W, I denoting a list of integers indices I = 1,2,3, •••, n of vessel trajectories in a set of vessel trajectories T, s denoting a number of vessel trajectories in a cluster, h denoting a vessel trajectory sum in a cluster which being a K × 2 matrix and being equal to Eq. (7): $\text{h}={\displaystyle {.Math.}_{\text{i=1}}^{\text{i=s}}{\text{T}}_{\text{i}}}$ wherein T.sub.i denoting a K × 2 matrix of an ith vessel trajectory, ${\displaystyle {.Math.}_{\text{i=1}}^{\text{i=s}}{\text{T}}_{\text{i}}}$ denoting a matrix summation; denoting a centroid vessel trajectory v as shown in Eq. (8): $\text{v}\text{=}\text{h}/\text{s}$ denoting a direct distance d.sub.d, a flip distance d.sub.F and a minimum average direct-flip distance MDF as shown in Expression set (9): $\left\{\begin{array}{c}{\text{d}}_{\text{d}}\left(\text{P,Q}\right)=\frac{1}{\text{k}}{\displaystyle \underset{\text{i}=1}{\overset{\text{k}}{.Math.}}\left|{\text{P}}_{\text{i}}-{\text{Q}}_{\text{i}}\right|}\\ {\text{d}}_{\text{F}}\left(\text{P,Q}\right)=\text{d}\left(\text{P,Q}\right)=\text{d}\left({\text{P}}_{\text{F}}\text{,Q}\right)\\ \text{MDF}\left(\text{P,Q}\right)=\min \left({\text{d}}_{\text{d}}\left(\text{P,D}\right),{\text{d}}_{\text{F}}\left(\text{P,Q}\right)\right)\end{array}\right)$ wherein |P.sub.i - Q.sub.i| denoting a distance between vessel trajectory point P.sub.i and vessel trajectory point Q.sub.i, the direct distance d.sub.d(P,Q) between two vessel trajectories denoting an mean distance between corresponding points of vessel trajectory P and vessel trajectory Q, a flip distance d.sub.F(P,Q) denoting a mean distance between a vessel trajectory and a corresponding points of another vessel trajectory after the flip, and the minimum average direct-flip distance MDF(P, Q) denoting a minimum of the direct distance d.sub.d(P,Q) and the flip distance d.sub.F(P,Q); initiating as follows: selecting a first vessel trajectory T.sub.1 and putting it to a first cluster c.sub.1, W = 1, C = {c.sub.1}, c.sub.1 = ({1}, T.sub.1, 1), obtaining a centroid vessel trajectory v.sub.1 = T.sub.1 in the first cluster c.sub.1 by Eq. (8), for each remaining vessel trajectories in turn T = {T.sub.i}, i = 2,3, ..., n which a total number of n — 1 vessel trajectories: calculating average direct-flip distances MDF(v.sub.1, T.sub.i) between remaining vessel trajectories T.sub.i and a centroid vessel trajectory v.sub.1 with Expression set (9), adding a vessel trajectory T.sub.d with a minimum value MDF(v.sub.1, T.sub.d) in MDF(v.sub.1, T.sub.i)to the first cluster c.sub.1 if any average minimum direct flip distances MDF(v.sub.1, T.sub.i) being less than a clustering threshold σ, obtaining c.sub.1 = ({1, d}, T.sub.1 + T.sub.d, 1 + 1) and ${\text{v}}_{\text{1}}=\frac{{\text{T}}_{\text{1}}{\text{+T}}_{\text{d}}}{2}$ in the first cluster c.sub.1, for each remaining vessel trajectories in turn T = {T.sub.i}, i = 2,3, ..., n which a total number of n — 2 vessel trajectories, processing each remaining vessel trajectories T.sub.i by step (4.3); otherwise creating a new cluster c.sub.2, selecting a vessel trajectory T.sub.d with a minimum value MDF(v.sub.1, T.sub.d) greater than the clustering threshold σ, c.sub.2 = ({d}, T.sub.d, 1), C = {c.sub.1, c.sub.2}, for each remaining vessel trajectories in turn T.sub.i = {T.sub.2, T.sub.3, ..., T.sub.n} which a total number of n — 2 vessel trajectories, processing each remaining vessel trajectories T.sub.i by step (4.3); (4.3) calculating minimum direct flip distances MDF(v.sub.e, T.sub.i) between remaining vessel trajectories T.sub.i and a centroid vessel trajectory v.sub.e of all the current clusters c.sub.e, e = 1, ... W with Expression set (9); adding vessel trajectory T.sub.i to a cluster c.sub.e with a minimum value for MDF(v.sub.e, T.sub.i), c.sub.e = ({I,i},h + T.sub.i, s + 1) if any average minimum direct flip distances MDF(v.sub.e, T.sub.i)being less than a clustering threshold σ; otherwise creating a new cluster c.sub.W+1, c.sub.W+1 = ({i}, T.sub.i, 1), incrementing W by 1; continuing to process steps (4.3) for remaining vessel trajectories T.sub.i in T until T={ }.

Description

BRIEF DESCRIPTION OF DRAWINGS

[0033] In order to illustrate the technical solution of the invention more clearly, the following is a brief description of the accompanying drawings to be used in the description, and it is obvious that the following drawings in the description are embodiments of the invention, from which other drawings can be obtained without creative work for a person of ordinary skill in the art.

[0034] FIG. 1 is schematic diagram of overall process of the method for vessel traffic pattern recognition via data quality control and data compression of the present invention.

[0035] FIG. 2 is a schematic diagram of a single vessel trajectory compression process of the method for vessel traffic pattern recognition via data quality control and data compression of the present invention.

[0036] FIG. 3 is a schematic diagram of Douglas-Peucker Pseud-Code process for a single vessel trajectory of the method for vessel traffic pattern recognition via data quality control and data compression of the present invention.

[0037] FIG. 4 is a schematic diagram of Quick Bundles algorithm clustering process of the method for vessel traffic pattern recognition via data quality control and data compression of the present invention.

[0038] FIG. 5 is a schematic diagram of Quick Bundles Pseud-Code process of the method for vessel traffic pattern recognition via data quality control and data compression of the present invention.

[0039] FIG. 6 is an original voyage trajectory of the method for vessel traffic pattern recognition via data quality control and data compression of the present invention.

[0040] FIG. 7 is a vessel’s repaired trajectory of the method for vessel traffic pattern recognition via data quality control and data compression of the present invention, with the dot in the figure showing the missing location of the trajectory detected and repaired based on the AIS update mechanism.

[0041] FIG. 8 is a total average compression rate and a total compression error under different compression thresholds of the method for vessel traffic pattern recognition via data quality control and data compression of the present invention.

[0042] FIG. 9 is a pre-compression vessel trajectory of the method for vessel traffic pattern recognition via data quality control and data compression of the present invention.

[0043] FIG. 10 is a compressed vessel trajectory of the method for vessel traffic pattern recognition via data quality control and data compression of the present invention.

[0044] FIG. 11 is a type of vessel trajectory similarity metric in the same direction of the method for vessel traffic pattern recognition via data quality control and data compression of the present invention.

[0045] FIG. 12 is a type of ship trajectory similarity metric in the reverse direction of the method for vessel traffic pattern recognition via data quality control and data compression of the present invention.

[0046] FIG. 13 shows major movement patterns of the vessel in the study area in step (4) of the preferred embodiment of the method for vessel traffic pattern recognition via data quality control and data compression of the present invention.

EMBODIMENTS

[0047] In order to better understand the technical features, objectives and effects of the present invention, the invention is described in more detail below in conjunction with the accompanying drawings. It should be understood that the specific embodiments described herein are intended to explain the invention only and are not intended to limit the patent of the invention. It should be noted that these drawings are in a very simplified form and use non-precise ratios only to facilitate and clearly assist in illustrating the patent of the invention.

[0048] A vessel traffic pattern recognition method incorporating data quality control and data compression is shown in FIG. 1 and includes the following steps:

[0049] assorting a collection of AIS data points according to MMSI and sorting each collection result by time ascending order to achieve stripping of AIS data points from different vessels: allocating each AIS data point in a collection to a vessel trajectory trajectory.sub.z so that each AIS data point therein having a same MMSI, sorting each vessel trajectory trajectory.sub.z by time ascending order, thus obtaining a set of vessel trajectories trajectory = {trajectory.sub.z}, z = 1,2,3, ...,243.

[0050] In the embodiment, each AIS data point of a vessel trajectory trajectory.sub.z represented by e = {MMSI, Time, lon, lat, sog}, MMSI denote a Maritime Mobile Service Identify of vessel, Time denote a time stamp, lon denote a longitude, lat denote a latitude, and sog denote a vessel speed over ground for said each vessel trajectory trajectory.sub.z.

[0051] A total of 243 vessel trajectories were collected and a partial information of trajectory.sub.1 is shown in Table 1.

TABLE-US-00001 partial information of trajectory.sub.1 MMSI Time lon lat sog 412358280 2019/1½ 7:35 122.2006 30.71712 8.4 412358280 2019/1½ 7:36 122.2006 30.716 8.2 412358280 2019/11/20 11:13 122.1433 30.52977 6.5 412358280 2019/11/20 11:14 122.1419 30.52839 6.6

Deleting duplicative AIS data points and segmenting vessel trajectory for each vessel trajectory trajectory.sub.z as following: for AIS data points therein having a same time stamp, a same longitude, a same latitude, and a same vessel speed over ground retaining only one thereof, while deleting the others thereof, thereafter segmenting vessel trajectory, starting from index 1 in trajectory.sub.z to obtain a first AIS data point efirst(j - 1) and a last AIS data point elast(j) such that AIS data points therebetween satisfying constraint in Expression set (1), continuing till end of index of trajectory.sub.z while deleting all the AIS data points between the first AIS data point efirst(j - 1) and the last AIS data point elast(j), and segmenting vessel trajectory trajectory.sub.z with elast(j) as a AIS data first point of a trajectory segment tra.sub.i, obtaining a new set of vessel trajectories tra = {tra.sub.i}, i = 1,2,3, ... 403 , wherein tra¡ denoting a i th vessel trajectory which i = 1,2,3, ... 403 , each AIS data point of a vessel trajectory tra¡ represented by e = {MMSI, Time, lon, lat, sog}.

[00012]$\left\{\begin{array}{c}{\text{sog}}_{\text{j}}<1\\ {\text{time}}_{\text{elast}\left(\text{j}\right)}-{\text{time}}_{\text{efirst}\left(\text{j}-1\right)}>{\text{Time}}_{\max }\end{array}\right)$

wherein sog.sub.j denoting a speed over ground at a jth AIS data point in a vessel trajectory, time.sub.efjrst(j-1) denoting a timestamp of an AIS data point efirst(j - 1) in a vessel trajectory, time.sub.elast(j) denoting a timestamp of an AIS data point elast(j) in a vessel trajectory, and Time.sub.max denoting a set time threshold.

[0052] In the embodiment, data from a total of 243 vessels are processed, after vessel trajectory segmentation process, 403 valid vessel trajectories are obtained.

[0053] Identifying adrift AIS data points and missing vessel trajectory segments for each vessel trajectory, repairing the missing vessel trajectory segments with cubic spline interpolation algorithm after deleting the adrift AIS data points to obtain high-quality AIS data, steps for each vessel trajectory tra¡ are as follows: [0054] (2.1) Setting a maximum safe driving speed of 30 knots, calculating a maximum displacement Δd.sub.j of adjacent AIS data points e.sub.j-1 to e.sub.j and a maximum displacement Δd.sub.j+1 of adjacent AIS data points e.sub.j to e.sub.j+1 according to the maximum safe driving speed of 30 knots to obtain a maximum longitude displacement value and a maximum latitude displacement value of adjacent AIS data points e.sub.j-1 to e.sub.j and e.sub.j to e.sub.j+1, calculating a longitude displacement difference Δlon.sub.j and a latitude displacement difference Δlat.sub.j from e.sub.j-1 to e.sub.j and a longitude displacement difference Δlon.sub.j+1 and a latitude displacement difference Δlat.sub.j+1 from e.sub.j to e.sub.j+1 respectively, a AIS point e.sub.j being a adrift AIS point if the longitude displacement difference Δlon.sub.j, Δlon.sub.j+1 and the latitude displacement difference Δlat.sub.j, Δlat.sub.j+1 satisfying a constraint of Expression set (2), and deleting the adrift AIS point e.sub.j; wherein Δt.sub.j denoting a time interval from adjacent AIS data points e.sub.j-1 to e.sub.j in a vessel trajectory, Time.sub.j-1 denoting a time stamp of an AIS data point e.sub.j-1, Time.sub.j denoting a time stamp of an AIS data point e.sub.j, Δt.sub.j+1 denoting a time interval from adjacent AIS data points e.sub.j+1 to e.sub.j in a vessel trajectory, Time.sub.j+1 denoting a time stamp of an AIS data point e.sub.j+1; [0055] (2.2) identifying missing vessel trajectory segments, a vessel trajectory of adjacent AIS data points will be regarded as a trajectory missing segment if a time interval between adjacent AIS data points is greater than 3 min but less than 5 min; [0056] (2.3) repairing the missing vessel trajectory segments by cubic spline interpolation algorithm in Eq. (4) subsequent to deletion of the adrift AIS data points in step (2.1) to obtain high-quality AIS data, for each missing vessel trajectory segment as follows: dividing a time series [A, B] of missing vessel trajectory segment into u intervals according to a time interval of 30 seconds, namely [[x.sub.1, x.sub.2], [x.sub.2, x.sub.3], ..., [x.sub.u, x.sub.u+1]], each sub-time series [x.sub.1, x.sub.2], [x.sub.2,x.sub.3], ..., [x.sub.u-1,x.sub.u] with 30 seconds time interval, a time interval of a sub-time series [x.sub.u, x.sub.u+1] being less than or equal to 30 seconds, A ≤ x.sub.1 < x.sub.2 < ••• < x.sub.u < x.sub.u+1 ≤ B; x.sub.1,x.sub.2,x.sub.3, ..., x.sub.u+1 corresponding to function values of y.sub.1,y.sub.2,y.sub.3, ...,y.sub.u+1 with y.sub.U = S(x.sub.U), (U = 1,2, ...,u), each sub-time series [x.sub.U, x.sub.U+1] satisfying Eq. (4); interpolating a longitude lon and a latitude lat and a vessel speed over ground sog of each time point x.sub.U in the missing vessel trajectory segment, y denoting a longitude lon when interpolating a longitude of a time point, y denoting a latitude lat when interpolating a latitude of a time point, y denoting a vessel speed over ground sog when interpolating a vessel speed over ground of a time point, obtaining a new vessel track.sub.i after a vessel trajectory repair; wherein a.sub.U, b.sub.U, c.sub.U, d.sub.U denoting pending coefficients which being derived from the missing vessel trajectory segment; In the embodiment, processing 403 vessel trajectories are processed to identify 3089 adrift AIS data points and 365 missing vessel trajectory segments, obtaining a new set of vessel trajectories track = {track.sub.i}, i = 1,2,3, ... 403, subsequent to processing of each vessel trajectory tra.sub.i in step (2), wherein track.sub.i denotes an ith vessel trajectory for i = 1,2,3, ... 403, each AIS data point of a vessel trajectory track.sub.i represented by e = {MMSI, Time, lon, lat, sog} after deleting the adrift AIS data points and repairing missing segments of vessel trajectory by cubic spline interpolation algorithm. The interpolation effect is shown in FIG. 6 and FIG. 7. The dots shown in FIG. 7 are interpolation points. The above effectively identifies and repairs the abnormal data in the vessel trajectory.

[0057] Compressing each vessel trajectory track.sub.i with a Douglas-Peucker algorithm by means of a self-invoking computer program as step (3.3) (reducing computational expenses in the clustering process of step (4)), as follows:

[0058] In the embodiment, to determine an optimal compression threshold of the Douglas-Peucker algorithm, testing a compression effect of the Douglas-Peucker algorithm under a compression threshold of 0 m, 0.5 m, ..., 20 m respectively, a compression rate being 71.4% and a compression error reaches 1.3 m when increasing the compression threshold to 12 m; with increasing the compression threshold further, a compression rate of the vessel trajectory data changes slowly, but a compression error of the data increases sharply; considering factors such as compression ratio and compression error, setting the compression threshold to 12 m in this embodiment, when the compression threshold being 12 m, a compression ratio being 44% and a compression error being 1.93 m. A total average compression ratio and total compression error under different compression thresholds is shown in FIG. 8. According to the compression threshold 12 m, a compression steps for each vessel trajectory track.sub.i are as follows:

[0059] (3.1) Forming a set of vessel trajectory points p = {p.sub.j(lon.sub.j,lat.sub.j)},j = 1,2,3, ..., v from the vessel trajectory track.sub.i, wherein p.sub.j denoting a jth vessel trajectory point for j = 1,2,3, ...,v, lon.sub.j denoting a jth longitude value in vessel trajectory point p.sub.j, lat.sub.j denoting a jth latitude value in vessel trajectory point p.sub.j; converting each vessel trajectory point p.sub.j from longitude and latitude coordinates to a Mercator coordinates vessel trajectory point m.sub.j with Equation set (5), thus obtaining M = {m.sub.j (mlon.sub.j, mlat.sub.j)},j = 1,2,3, ...,v, wherein M denoting a set of vessel trajectory points in the Mercator coordinate system and M = {m.sub.1(mlon.sub.1, mlat.sub.1), m.sub.2(mlon.sub.2, mlat.sub.2), m.sub.3(mlon.sub.3, mlat.sub.3), ..., m.sub.v(mlon.sub.v, mlat.sub.v)} , m.sub.j denoting a jth vessel trajectory point in the Mercator coordinate system which j = 1,2,3, ...,v, mlon.sub.j denoting a jth longitude value in vessel trajectory point m.sub.j in Mercator coordinate system, mlat.sub.j denoting a jth latitude value in vessel trajectory point m.sub.j in the Mercator coordinate system;

[00016]$\left\{\begin{array}{c}\text{radius}=\frac{\text{lr}\ast \cos \beta }{\sqrt{1-{\text{E}}^2\ast {\sin }^2\beta }}\\ {\text{q}}_{\text{j}}=\ln \left(\tan \left(\frac{\text{π}}{4}+\frac{{\text{lat}}_{\text{j}}}{2}\right){\left(\frac{1-\text{E}\ast \sin {\text{lat}}_{\text{j}}}{1+\text{E}\ast \sin {\text{lat}}_{\text{j}}}\right)}^2\right)\\ {\text{Mlon}}_{\text{j}}=\text{radius}\ast {\text{lon}}_{\text{j}}\\ {\text{Mlat}}_{\text{j}}=\text{radius}\ast {\text{q}}_{\text{j}}\end{array}\right)$

[0060] wherein radius denoting a radius of the standard latitude-parallel circle, lr denoting a long radius of Earth’s ellipsoid, β a standard latitude in the Mercator projection, E denoting a first eccentricity of Earth’s ellipsoid, q.sub.j denoting an equivalent latitude of a jth vessel trajectory point; (3.2) initiating in respective of the set of vessel trajectory points M = {m.sub.1(mlon.sub.1, mlat.sub.1), m.sub.2(mlon.sub.2, mlat.sub.2), m.sub.3(mlon.sub.3, mlat.sub.3), ..., m.sub.v(mlon.sub.v, mlat.sub.v)} as follows: denoting r as a set of key vessel trajectory points, putting a starting vessel trajectory point m.sub.1(mlon.sub.1, mlat.sub.1) and an end vessel trajectory point m.sub.v(mlon.sub.v, mlat.sub.v) in the set of vessel trajectory points M as key vessel trajectory points to the set of key vessel trajectory points r in order, obtaining r = {m.sub.1(mlon.sub.1, mlat.sub.1, m.sub.v(mlon.sub.v, mlat.sub.v)}; connecting the starting vessel trajectory point m.sub.1(mlon.sub.1, mlat.sub.1) and the end vessel trajectory point m.sub.v(mlon.sub.v, mlat.sub.v) in the set of vessel trajectory points M as a straight line l.sub.1v , calculating distances dist = {dist.sub.2, dist.sub.3, ..., dist.sub.v-1} from all vessel trajectory points between m.sub.1(mlon.sub.1, mlat.sub.1) and m.sub.v(mlon.sub.v, mlat.sub.v) to the straight line l.sub.1v, with Eq. (6), determining a vessel trajectory point m.sub.g(mlon.sub.g, mlat.sub.g) such that dist.sub.g = max {dist.sub.2, dist.sub.3, ..., dist.sub.v-1}; [0061] wherein dist denoting a vertical distance from a vessel trajectory point to a straight line in the Mercator coordinate system, se denoting a vector from a start of the straight line to an end of the straight line, ta denoting a vector from the start of the straight line to a target point; [0062] wherein dist denoting a vertical distance from a vessel trajectory point to a straight line in the Mercator coordinate system, se denoting a vector from a start of the straight line to an end of the straight line, ta denoting a vector from the start of the straight line to a target point; [0063] concluding step(3.2) on condition dist.sub.g being less than a set compression threshold 12 m; otherwise, putting the vessel trajectory point m.sub.g(mlon.sub.g, mlat.sub.g) as a key vessel trajectory point to r in order, obtaining r = {m.sub.1(mlon.sub.1, mlat.sub.1), m.sub.g(mlon.sub.g, mlat.sub.g), m.sub.v(mlon.sub.v, mlat.sub.v)} , dividing the set of vessel trajectory points M = {m.sub.1(mlon.sub.1, mlat.sub.1), m.sub.2(mlon.sub.2, mlat.sub.2), m.sub.3(mlon.sub.3, mlat.sub.3), ..., m.sub.v(mlon.sub.v, mlat.sub.v)} into two sub vessel trajectory point sets M.sub.gsub.sub.h,h = 1,2 from m.sub.1(mlon.sub.1, mlat.sub.1) to m.sub.g(mlon.sub.g, mlat.sub.g) and m.sub.g(mlon.sub.g, mlat.sub.g) to m.sub.v(mlon.sub.v, mlat.sub.v) , M.sub.gsub.sub.1 = {m.sub.1(mlon.sub.1, mlat.sub.1), ..., m.sub.g(mlon.sub.g, mlat.sub.g)} from m.sub.1(mlon.sub.1, mlat.sub.1) to m.sub.g(mlon.sub.g, mlat.sub.g) and M.sub.gsub.sub.2 = {m.sub.g(mlon.sub.g, mlat.sub.g), ..., m.sub.v(mlon.sub.v, mlat.sub.v)} form m.sub.g(mlon.sub.g, mlat.sub.g) to m.sub.v(mlon.sub.v, mlat.sub.v), wherein M.sub.gsub.sub.1 denoting a first set of sub vessel trajectory points, M.sub.gsub.sub.2 denoting a 2nd set of sub vessel trajectory points; calculating a number of vessel trajectory points M.sub.gsub.sub.1number.sub.1 in M.sub.gsub.sub.1 and a number of vessel trajectory points M.sub.gsub.sub.1number.sub.2 in M.sub.gsub.sub.2 , processing M.sub.gsub.sub.1 by step (3.3) if the number of vessel trajectory points M.sub.gsub.sub.1number.sub.1 being greater than a set number threshold 50; processing M.sub.gsub.sub.2 by step (3.3) if the number of vessel trajectory points M.sub.gsub.sub.1number.sub.2 being greater than the set number threshold 50;

[0064] (3.3) Mtrack = {m.sub.start(mlon.sub.start, mlat.sub.start), ..., m.sub.end(mlon.sub.end, mlat.sub.end)} denoting a sub vessel trajectory point set, m.sub.start(mlon.sub.start, mlat.sub.start) denoting a first vessel trajectory point which start = 1,2,3, ...,v - 1, m.sub.end(mlon.sub.end, mlat.sub.end) denoting a last vessel trajectory point which end = 2,3, ..., v, a subscript start being less than subscript point end; connecting the first point m.sub.start(mlon.sub.start, mlat.sub.start) and the last point m.sub.end(mlon.sub.end, mlat.sub.end) as a straight line l.sub.startend, calculating distances dist = {dist.sub.start+1, dist.sub.start+2, ..., dist.sub.end-1,} from all vessel trajectory points between m.sub.start(mlon.sub.start, mlat.sub.start) and m.sub.end(mlon.sub.end, mlat.sub.end) to the straight line l.sub.startend with Eq. (6), determining a vessel trajectory point m.sub.d(mlon.sub.d, mlat.sub.d) such that dist.sub.d = max{dist.sub.start+1, dist.sub.start+2, ..., dist.sub.end-1}, concluding step (3.3) on condition dist.sub.d being less than the compression threshold 12 m; otherwise, putting the vessel trajectory point m.sub.d(mlon.sub.d, mlat.sub.d) as a key vessel trajectory point to r, dividing the sub vessel trajectory point set Mtrack into two sub vessel trajectory point sets M.sub.dsub.sub.h,h = 1,2 from m.sub.start(mlon.sub.start, mlat.sub.start) to m.sub.d(mlon.sub.d, mlat.sub.d) and m.sub.d(mlon.sub.d, mlat.sub.d) to m.sub.end(mlon.sub.end, mlat.sub.end), M.sub.dsub.sub.1 = {m.sub.start(mlon.sub.start, mlat.sub.start), ..., m.sub.d(mlon.sub.d, mlat.sub.d)} and M.sub.dsub.sub.2 = {m.sub.d(mlon.sub.d, mlat.sub.d), ..., m.sub.end(mlon.sub.end, mlat.sub.end)}, wherein M.sub.dsub.sub.1 denoting a first set of sub vessel trajectory points after splitting the sub vessel trajectory point set Mtrack with the vessel trajectory point m.sub.d(mlon.sub.d, mlat.sub.d) as a split point, M.sub.dsub.sub.2 denoting a 2nd set of sub vessel trajectory points after splitting the sub vessel trajectory point set Mtrack with the vessel trajectory point m.sub.d(mlon.sub.d, mlat.sub.d) as a split point; calculating a number of vessel trajectory points M.sub.dsub.sub.1number.sub.1 in M.sub.dsub.sub.1 and a number of vessel trajectory points M.sub.dsub.sub.1number.sub.2 in M.sub.dsub.sub.2 , processing M.sub.dsub.sub.1 by step (3.3) if the number of vessel trajectory points M.sub.dsub.sub.1number.sub.1 being greater than a set number threshold 50, processing M.sub.dsub.sub.2 by step (3.3) if the number of vessel trajectory points M.sub.dsub.sub.1number.sub.2 being greater than the set number threshold 50 until the subscript start greater being than or equal to end.

[0065] In the embodiment, processing 403 vessel trajectories to obtain a new set of vessel trajectories R = {r.sub.i}, i = 1,2,3, ... 403, wherein r.sub.i denoting a vessel trajectory of ith vessel which i = 1,2,3, ... 403 , each vessel trajectory points of vessel trajectory r.sub.i represented by m = {mlon, mlat}. A schematic diagram of a single vessel trajectory compression process is shown FIG. 2. Douglas-Peucker Pseudo-Code for a vessel trajectory is shown in Table 2. A schematic diagram of Douglas-Peucker Pseud-Code process for a single vessel trajectory is shown FIG. 3. The effect of a single vessel voyage trajectory before compression is shown in FIG. 9, and the effect after compression is shown in FIG. 10.

TABLE-US-00002 Douglas-Peucker Pseudo-Code for a vessel trajectory Algorithm: Douglas-Peucker Pseudo-Code Input: a set of trajectory points of a vessel trajectory m = {m.sub.1,m.sub.2,m.sub.3, ..., m.sub.v} 1:index = 1 2: end = len(m) 3. def compression (self, m, start, endpoint): 4: r= {m.sub.1, m.sub.v} # r denotes a set of key vessel trajectory points 5: if len(m[start: endpoint]) > .Math. then # .Math. denotes a set number threshold 6: d.sub.max = 0 :7 currentIndex = 1 8: for i in range(start + 1, endpoint - 1) do 9: distance = dist(m.sub.i, line(m.sub.start,m.sub.endpoint)) 10 if distance > d.sub.max then 11: d.sub.max = distance 12: currentIndex = i 13: if d.sub.max > ε then # ε denotes a set compression threshold 14: append (r, m.sub.i) 15: self. compression (m, start, currentIndex) 16: self. compression (m, currentIndex, endpoint) 17: return r 18: r = compression (m, index, end) Output: r

[0066] Reconstructing each vessel trajectory r.sub.i with cubic spline interpolation algorithm, and clustering vessel trajectories into various clusters by Quick Bundles algorithm to form a vessel traffic pattern as follows:

[0067] (4.1) reconstructing each vessel trajectory r.sub.i with cubic spline interpolation algorithm, for each vessel trajectory r.sub.i in R, searching a vessel trajectory r.sub.j with most vessel trajectory points, calculating number differences between remaining vessel trajectories and the vessel trajectory r.sub.j trajectory points respectively, and interpolating at the end of each remaining vessel trajectory with cubic spline interpolation algorithm so that each vessel trajectory has same number of trajectory points to obtain a new set of vessel trajectories T = {T.sub.i{t.sub.j(mlon.sub.j, mlat.sub.j)|j = 1,2,3, ... ,4578}],i = 1,2,3, ... 403, wherein T.sub.i denoting an ith vessel trajectory which i = 1,2,3, ... 403, each vessel trajectory T.sub.i being a 4578 × 2 matrix; t.sub.j denoting an jth vessel trajectory point of time order serial number j = 1,2,3, ...,4578, each vessel trajectory point t.sub.j of a vessel trajectory T.sub.i represented by t = {mlon, mlat}; each vessel trajectory T.sub.i = (t.sub.1,t.sub.2, .Math..Math..Math., t.sub.4578) has two ordered polylines, namely a isotropic trajectory T.sub.i = (t.sub.1, t.sub.2, .Math..Math..Math. t.sub.4578) and a reverse trajectory flip version T.sub.Fi = (t.sub.4578, t.sub.4578-1, .Math..Math..Math. t.sub.1);

[0068] (4.2) clustering vessel trajectories into various clusters by Quick Bundles algorithm to form a vessel traffic pattern: constructing a cluster class set of vessel trajectories C = {c.sub.q(I, h, s)|q = 1,2, ..., W}, wherein c.sub.q denoting a cluster set of vessel trajectories in cluster q which q = 1,2, ..., W, I denoting a list of integers indices I = 1,2,3, ... ,403 of vessel trajectories in a set of vessel trajectories T, s denoting a number of vessel trajectories in a cluster, h denoting a vessel trajectory sum which being a 4578 × 2 matrix and being equal to Eq. (7):

[00018]$\text{h}={\displaystyle {.Math.}_{\text{i}=1}^{\text{i}=\text{s}}{\text{T}}_{\text{i}}}$

[0069] wherein T.sub.i denoting a 4578 × 2 matrix of an ith vessel trajectory, denoting a matrix summation; [0070] denoting a centroid vessel trajectory v as shown in Eq. (8): [0071] denoting a direct distance d.sub.d, a flip distance d.sub.F and a minimum average direct-flip distance MDF as shown in Expression set (9): [0072] wherein |P.sub.i - Q.sub.i| denoting a distance between vessel trajectory point P.sub.i and vessel trajectory point Q.sub.i, a direct distance d.sub.d(P, Q) between two trajectories denoting an mean distance between corresponding points of vessel trajectory P and vessel trajectory Q, a flip distance d.sub.F(P,Q) denoting a mean distance between a vessel trajectory and a corresponding points of another vessel trajectory after the flip, and a minimum direct flip distance MDF(P, Q) denoting a minimum of the direct distance d.sub.d(P,Q) and the flip distance d.sub.F(P,Q); [0073] In the embodiment, calculating a similarity matrix between vessel trajectories uses Equation set (9), a schematic diagram of vessel trajectory similarity metric type is shown in FIG. 11 and FIG. 12. Initiating as follows: selecting a first vessel trajectory T.sub.1 and putting it to a first cluster c.sub.1, W = 1, C = {c.sub.1}, c.sub.1 = ({1}, T.sub.1, 1), obtaining a centroid vessel trajectory v.sub.1 = T.sub.1 in the first cluster c.sub.1 by Eq. (8), for each remaining vessel trajectories in turn T = {T.sub.i},i = 2,3, ...,403 which a total number of 402 vessel trajectories: calculating minimum direct flip distances MDF(v.sub.1, T.sub.i) between remaining vessel trajectories T.sub.i and a centroid vessel trajectory v.sub.1 with Equation set (9), adding a vessel trajectory T.sub.d with a minimum value MDF(v.sub.1, T.sub.d) in MDF(v.sub.1, T.sub.i) to the first cluster c.sub.1 if any minimum direct flip distances MDF(v.sub.1, T.sub.i) being less than a clustering threshold σ, obtaining c.sub.1 = ({1, d}, T.sub.1 + T.sub.d, 1 + 1) and in the first cluster c.sub.1, number of remaining vessel trajectories being 401, processing each remaining vessel trajectories T.sub.i by step (4.3); otherwise creating a new cluster c.sub.2, selecting the vessel trajectory T.sub.d with a minimum value MDF(v.sub.1, T.sub.d) greater than the clustering threshold σ, c.sub.2 = ({d}, T.sub.d, 1), C = {c.sub.1, c.sub.2}, number of remaining vessel trajectories being 401, processing each remaining vessel trajectories T.sub.i by step (4.3);

[0074] (4.3) calculating minimum direct flip distances MDF(v.sub.e, T.sub.i) between remaining vessel trajectories T.sub.i and a centroid vessel trajectory v.sub.e of all the current clusters c.sub.e, e = 1, ... M with Equation set (9); adding vessel trajectory T.sub.i to a cluster c.sub.e with a minimum value for MDF(v.sub.e, T.sub.i), c.sub.e = ({I,i},h + T.sub.i, s + 1) if any minimum direct flip distances MDF(v.sub.e, T.sub.i) being less than a clustering threshold σ; otherwise creating a new cluster c.sub.M+1, c.sub.M+1 = ({i}, T.sub.i,1), incrementing M by 1, continuing to process steps (4.3) for remaining vessel trajectories T.sub.i in T until T={ }.

[0075] In the embodiment, 403 vessel trajectories are processed, and are clustered into various clusters by Quick Bundles algorithm to form vessel traffic patterns. A schematic diagram of Quick Bundles algorithm clustering process is shown in FIG. 4. A pseudo-code for Quick Bundles algorithm is shown in Table 4. A schematic diagram of Quick Bundles Pseud-Code process is shown in FIG. 5. A resulting cluster is shown in Table 3. A visualization effect of clustering of this implementation is shown in FIG. 13.

[0076] The dataset utilized therefor was collected in Shanghai Yangshan Port in a rectangle from (121.94◦E, 30.52◦N) to (122.22◦E, 30.72◦N) were analyzed, comprising AIS observations of vessels from Nov. 01, 2019 to Nov. 30, 2019. The raw dataset contains 1,004,121 pieces of AIS data points. The patterns displayed in FIG. 13 show that: a majority of vessels are more active in the southwest side of Xiaoyangshan deep-water port area and an east side of Bojiazui Island, while relatively few vessels are in the north side of Xiaoyangshan deep-water port area or the northeast side of Little Turtle Island. The results of the embodiment prove the feasibility of the present invention in understanding vessel traffic patterns for maritime factual real-time supervision and in discovering distribution of vessel trajectory activities among scattered and chaotic vessel traffic.

[0077] As can be seen thereabove, steps (1), (2), and (3) are pre-processing steps for processing the raw AIS data, that is, the collection of AIS data points, to obtain a set of vessel trajectories as below: T = {T.sub.j{t.sub.j(longitude.sub.j, latitude.sub.j)|j = k}}, wherein T.sub.i denote an ith vessel trajectory which i = 1,2,3, ... n, each vessel trajectory T.sub.i is a k × 2 matrix; t.sub.j denote an jth vessel trajectory point of time order serial number j = 1,2,3, ... k, each vessel trajectory point t.sub.j of a vessel trajectory T.sub.i represented by t = {longitude, latitude}. Thereafter, the afore-mentioned set of vessel trajectories is inputted into step (4) to obtain identification of the vessel traffic patterns. To conclude, step (4) per se works as an independent vessel trajectory clustering process for identification of the vessel traffic patterns.

TABLE-US-00003 Information of some vessel track segments after clustering Cluster category W MMSI Cluster class 1 219034000,219231000,.Math..Math..Math.,636017686,636018059 Cluster class 2 412254253,412371217,.Math..Math..Math.,412380360,413595000 Cluster class 3 412355690,412373080,.Math..Math..Math.,413304000,413557430 Cluster class 4 412358240,412358280,.Math..Math..Math.,413364330,413368640 Cluster class 5 412373080,412421040,.Math..Math..Math.,412373080,413557430

TABLE-US-00004 Quick Bundles Pseudo-Code Algorithm: Quick Bundles Pseudo-Code Input: T = {T.sub.1, T.sub.2, T.sub.3, ..., T.sub.n} 1: c.sub.1 = ([1], T.sub.1, 1) #creating first cluster 2: C = {c.sub.1} 3: W = 1 4: for i = 2 to n do 5: t=T.sub.i 6: alld=infinity(W) 7: flip=zeros(W) 8: for e=1 to W do 9: v = c.sub.eh/c.sub.es 10: d = d.sub.d(t,v) 11: f = d.sub.f(t,v) 12: if f<d then 13: d=f 14: flip=1 15: end if 16: alld=d 17: end for 18: m=min(alld) 19: 1=argmin(alld) 20: if m < σ then #σ denote a clustering threshold 21: if flip.sub.1 is 1 then 22: c.sub.1h = c.sub.1h + t.sub.f 23: else 24: c.sub.1h = c.sub.1h + t 25: end if 26: c.sub.1s = c.sub.1s + 1 27: append(c.sub.1I, i) 28: else 29: c.sub.W+1 = ([i], t, 1) 30: append(C, c.sub.W+1) 31: W=W+1 32: end if 33: end for Output: C = {c.sub.1, c.sub.2, c.sub.3, ..., c.sub.W}

[0078] As described above, it is only a specific embodiment of the present invention, but the scope of protection of the present invention is not limited to it, and any person skilled in the art can easily think of various equivalent modifications or substitutions within the scope of the technology disclosed herein, which shall be included in the scope of protection of the present invention. Therefore, the scope of protection of the present invention shall be subject to the scope of protection of the claims.