From e9f65e2c30a806adfd1bdff3b341b0730c49f49a Mon Sep 17 00:00:00 2001
From: amitjans <amitjans@ethz.ch>
Date: Mon, 31 May 2021 20:47:44 +0200
Subject: [PATCH] Compute 2d distance instead of 1d

---
 tools.py | 71 ++++++++++++++++++++++++++++++++------------------------
 1 file changed, 41 insertions(+), 30 deletions(-)

diff --git a/tools.py b/tools.py
index 9d79448..fe274d5 100644
--- a/tools.py
+++ b/tools.py
@@ -218,6 +218,14 @@ def ion_positions(image):
     Finally return the coordinates together with the maximum ion width.
     """
     coordinates, max_width = _max_coord(image, width=True) # measure local max
+    # delete local maxima that are not in the same axis as the brightest spot.
+    try:
+        h_axis = coordinates[0, 0]
+        coordinates = coordinates[np.logical_and(coordinates[:, 0] < h_axis + 15,
+        coordinates[:, 0] > h_axis - 15)]
+    except IndexError:
+        pass
+
     coordinates += np.flip(settings.REGION[0]) # change reference
     coordinates = coordinates[coordinates[:, 1].argsort()] # order
     # add 3rd column indicating the ion's state (bright, dark, wrong isotope)
@@ -287,8 +295,11 @@ def compensate_gradient(image, coord):
 
     slope = _measure_gradient(image, coord)
 
-    for i in range(image.shape[1]):
-        image[:, i] = np.clip(image[:, i] - i*slope, 0, None)
+    start = coord[0, 1] - 30
+    end = coord[-1, 1] + 30
+
+    for i in np.r_[start:end]:
+        image[:, i] = np.clip(image[:, i] - (i-start)*slope, 0, None)
 
     return image.astype(np.uint16)
 
@@ -335,7 +346,7 @@ def dark_ions(image, coord):
     """
     # Compute the positions starting from the left/right
     x_pos = np.stack((coord[:, 1] - coord[0, 1], np.sort(np.abs(coord[:, 1] - coord[-1, 1]))))
-    x_dist = _ion_distance(coord[:, 1])  # compute relative distances
+    x_dist = _ion_distance(coord)  # compute relative distances
     # Load the masks. TODO: Remove the mask options that end up outside the image!
     low_idx = max(0, len(x_pos[0]) - 2)
     high_idx = min(28, len(x_pos[0]) * 2)
@@ -364,33 +375,32 @@ def dark_ions(image, coord):
     comp[2, :, -1] = 0
     loss = np.sum(comp, axis=2)/len(x_pos)
 
-    if loss.min() <= 20:
-        idx = np.where(loss == np.min(loss))
-
-        if idx[0][0] == 1:  # if using "inverted positions", transform back
-            new_x_pos = np.abs(mask[idx[1][0], :][mask[idx[1][0], :] >= 0] - coord[-1, 1])
-        else:
-            new_x_pos = mask[idx[1][0], :][mask[idx[1][0], :] >= 0] + coord[0, 1]
-
-        diff_n = len(new_x_pos) - len(x_pos[0])
-        for i in range(diff_n):
-            mask2 = comparison[idx[0][0], idx[1][0], :] < 100000
-            idx2 = np.argsort(comparison[idx[0][0], idx[1][0], :][mask2])[-(diff_n - i)]
-            if idx2 != 0:
-                # If we are using the distances' loss, then we have to add 1 to the index
-                if idx[0][0] == 2:
-                    idx2 += 1
-                coord = np.append(coord, np.array([[np.mean(coord[:, 0]),
-                new_x_pos[idx2], 2]]), axis=0)
-
-            elif idx2 == 0 and idx[0][0] == 2:
-                coord = np.append(coord, np.array([[np.mean(coord[:, 0]), new_x_pos[0] -
-                                                        comparison[2, idx[1][0], 0], 2]]), axis=0)
-
-        coord = coord[coord[:, 1].argsort()].astype(int)
-        coord = coord[coord[:, 1] < image.shape[1]]  # remove points outside the image
+    idx = np.where(loss == np.min(loss))
 
+    if idx[0][0] == 1:  # if using "inverted positions", transform back
+        new_x_pos = np.abs(mask[idx[1][0], :][mask[idx[1][0], :] >= 0] - coord[-1, 1])
     else:
+        new_x_pos = mask[idx[1][0], :][mask[idx[1][0], :] >= 0] + coord[0, 1]
+
+    diff_n = len(new_x_pos) - len(x_pos[0])
+    for i in range(diff_n):
+        mask2 = comparison[idx[0][0], idx[1][0], :] < 100000
+        idx2 = np.argsort(comparison[idx[0][0], idx[1][0], :][mask2])[-(diff_n - i)]
+        if idx2 != 0:
+            # If we are using the distances' loss, then we have to add 1 to the index
+            if idx[0][0] == 2:
+                idx2 += 1
+            coord = np.append(coord, np.array([[np.mean(coord[:, 0]),
+            new_x_pos[idx2], 2]]), axis=0)
+
+        elif idx2 == 0 and idx[0][0] == 2:
+            coord = np.append(coord, np.array([[np.mean(coord[:, 0]), new_x_pos[0] -
+                                                    comparison[2, idx[1][0], 0], 2]]), axis=0)
+
+    coord = coord[coord[:, 1].argsort()].astype(int)
+    coord = coord[coord[:, 1] < image.shape[1]]  # remove points outside the image
+
+    if loss.min() > 20:
         warnings.warn("Loss is too high... "
         "Try calculating a new factor value using calibrate() function.")
 
@@ -543,7 +553,7 @@ def _max_coord(image, width=False):
     image = image[region[0, 1]:region[1, 1], region[0, 0]:region[1, 0]]
 
     filt_image = _difference_of_gaussians(image, low_sigma=settings.BLOB_SIZE, preserve_range=True)
-    coordinates = peak_local_max(filt_image, threshold_abs=0.85, min_distance=25)
+    coordinates = peak_local_max(filt_image, threshold_abs=0.85, min_distance=settings.ION_DIST)
 
     if width:
         try:
@@ -648,7 +658,8 @@ def _ion_distance(coord):
     """
     dist = []
     for i in range(len(coord) - 1):
-        dist.append(coord[i + 1] - coord[i])
+        dist.append(np.sqrt((coord[i, 0] - coord[i + 1, 0])**2 +
+        (coord[i, 1] - coord[i + 1, 1])**2).item())
     dist = np.array(dist)
 
     return dist
-- 
GitLab