diff --git a/skais/ais/ais_trajectory.py b/skais/ais/ais_trajectory.py
index dca3bcc230e78ae872fbcee55bcccc8709a54bdf..90ba3e83404f3ae828796b58d68027068f2cc538 100644
--- a/skais/ais/ais_trajectory.py
+++ b/skais/ais/ais_trajectory.py
@@ -1,5 +1,6 @@
 import numbers
 import random
+from copy import copy
 
 import pandas as pd
 import numpy as np
@@ -130,22 +131,9 @@ def generate_points(data, positions, height, width, node_size, lower_lat, upper_
                 data[x, y] = [1]
 
 
-@jit(nopython=True)
-def generate_links(data, positions, height, width, lower_lat, upper_lat, lower_lon, upper_lon):
-    lon, lat = positions[0, 0], positions[0, 1]
-    for longitude, latitude in positions[1:]:
-        x_prv, y_prev = get_coord(lat, lon, height, width, lower_lat, upper_lat, lower_lon, upper_lon)
-        x_nxt, y_nxt = get_coord(latitude, longitude, height, width, lower_lat, upper_lat, lower_lon,
-                                 upper_lon)
-
-        lon, lat = longitude, latitude
-        for x, y in bresenham(x_prv, y_prev, x_nxt, y_nxt):
-            data[x, y] = [1]
-
-
 @jit(nopython=True)
 def generate_points_with_features(data, positions, features_vectors, bounds, node_size, height, width,
-                                  lower_lat, upper_lat, lower_lon, upper_lon, link):
+                                  lower_lat, upper_lat, lower_lon, upper_lon):
     for pos, f in zip(positions, features_vectors):
         latitude = pos[1]
         longitude = pos[0]
@@ -162,21 +150,36 @@ def generate_points_with_features(data, positions, features_vectors, bounds, nod
             for y in range(y_lower_bound, y_upper_bound + 1):
                 for i, v in enumerate(value):
                     data[x, y, i] = v
-    if link:
-        lon, lat = positions[0, 0], positions[0, 1]
-
-        value = __get_image_value__(features_vectors[0], bounds)
-        for pos, f in zip(positions[1:], features_vectors[1:]):
-            latitude = pos[1]
-            longitude = pos[0]
-            x_prv, y_prev = get_coord(lat, lon, height, width, lower_lat, upper_lat, lower_lon, upper_lon)
-            x_nxt, y_nxt = get_coord(latitude, longitude, height, width, lower_lat, upper_lat, lower_lon,
-                                     upper_lon)
-            lon, lat = longitude, latitude
-            for x, y in bresenham(x_prv, y_prev, x_nxt, y_nxt):
-                for i, v in enumerate(value):
-                    data[x, y, i] = v
-            value = __get_image_value__(f, bounds)
+
+
+@jit(nopython=True)
+def generate_links(data, positions, height, width, lower_lat, upper_lat, lower_lon, upper_lon, values):
+    lon, lat = positions[0, 0], positions[0, 1]
+
+    current_value = values[0]
+    for pos, nxt_value in zip(positions[1:], values[1:]):
+        latitude = pos[1]
+        longitude = pos[0]
+        x_prv, y_prev = get_coord(lat, lon, height, width, lower_lat, upper_lat, lower_lon, upper_lon)
+        x_nxt, y_nxt = get_coord(latitude, longitude, height, width, lower_lat, upper_lat, lower_lon,
+                                 upper_lon)
+        lon, lat = longitude, latitude
+        for x, y in bresenham(x_prv, y_prev, x_nxt, y_nxt):
+            for i, v in enumerate(current_value):
+                data[x, y, i] = v
+        current_value = nxt_value
+
+
+@jit(nopython=True)
+def generate_values(features_vectors, bounds):
+    result = np.zeros(features_vectors.shape)
+    for i in range(len(features_vectors)):
+        features = features_vectors[i]
+        value = __get_image_value__(features, bounds)
+        for j in range(len(value)):
+            v = value[j]
+            result[i, j] = v
+    return result
 
 
 class AISTrajectory(AISPoints):
@@ -380,9 +383,14 @@ class AISTrajectory(AISPoints):
             generate_points(data, positions, height, width, node_size, lower_lat, upper_lat, lower_lon, upper_lon)
 
             if link:
-                generate_links(data, positions, height, width, lower_lat, upper_lat, lower_lon, upper_lon)
+                generate_links(data, positions, height, width, lower_lat, upper_lat, lower_lon, upper_lon,
+                               np.ones(len(positions)))
 
         else:
+            print(features_vectors[0] / 10)
             generate_points_with_features(data, positions, features_vectors, np.array(bounds), node_size, height, width,
-                                          lower_lat, upper_lat, lower_lon, upper_lon, link)
+                                          lower_lat, upper_lat, lower_lon, upper_lon, )
+            if link:
+                generate_links(data, positions, height, width, lower_lat, upper_lat, lower_lon, upper_lon,
+                               generate_values(features_vectors, bounds))
         return data
diff --git a/skais/tests/ais/test_ais_trajectory.py b/skais/tests/ais/test_ais_trajectory.py
index 51d776ccfe534ef89e1784ad5936dd68b8695572..1abe4b744ab5e401d3935676fab3fe88d4f67dc3 100644
--- a/skais/tests/ais/test_ais_trajectory.py
+++ b/skais/tests/ais/test_ais_trajectory.py
@@ -651,4 +651,31 @@ class TestAISTrajectoryImageGeneration(unittest.TestCase):
                              [[0.25, 1], [0, 0], [0, 0], [0, 0], [0, 0], [0, 0], [0, 0], [0, 0], [0, 0], [0, 0], [0, 0],
                               [0, 0], [0, 0], [0, 0], [0, 0], [0, 0], [0.25, 0.5], [0.5, 0.25]]])
 
-        np.testing.assert_array_equal(result, expected)
\ No newline at end of file
+        np.testing.assert_array_equal(result, expected)
+
+
+    def test_generate_array_interpolate_links(self):
+        trajectory = AISTrajectory(
+            pd.DataFrame(
+                {
+                    "latitude": [0, 10, 0, 20],
+                    "longitude": [0, 10, 20, 20],
+                    "ts_sec": [i for i in range(4)],
+                    "sog": [10, 10, 20, 40]
+                }
+            )
+        )
+
+        result = trajectory.generate_array_from_positions(height=9, width=18, link=True, bounding_box='fit',
+                                                          features={"sog": (0, 80)}, node_size=0).reshape((9, 18), interpolation='linear')
+        expected = np.array([[0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0.5],
+                             [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0.5],
+                             [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0.5],
+                             [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0.5],
+                             [0, 0, 0, 0, 0, 0, 0, 0.25, 0.25, 0.25, 0, 0, 0, 0, 0, 0, 0, 0.5],
+                             [0, 0, 0, 0, 0, 0.25, 0.25, 0, 0, 0, 0.25, 0.25, 0, 0, 0, 0, 0, 0.5],
+                             [0, 0, 0, 0.25, 0.25, 0, 0, 0, 0, 0, 0, 0, 0.25, 0.25, 0, 0, 0, 0.5],
+                             [0, 0.25, 0.25, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0.25, 0.25, 0, 0.5],
+                             [0.25, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0.25, 0.5]]) / 2
+
+        np.testing.assert_array_almost_equal(result, expected)