quatrope / astroalign

A tool to align astronomical images based on asterism matching
MIT License
140 stars 44 forks source link

Different results for parallel and sequential implementations of find_transform #85

Open prajwel opened 1 year ago

prajwel commented 1 year ago

Hi,

I noticed that I am getting different results for parallel and sequential implementations of find_transform. Please see the code below that reproduces the problem.

import numpy as np
import astroalign as aa

from joblib import Parallel, delayed

aa.NUM_NEAREST_NEIGHBORS = 6 

def get_rigid_transform(frame_sources, reference):
    transf, matched_positions = aa.find_transform(frame_sources, reference, max_control_points= 30)
    return transf.params

def match_frames_joblib_loop(sources_all_frames):
    trans_matrices = Parallel(n_jobs=-1)(delayed(get_rigid_transform)(i, sources_all_frames[0]) for i in sources_all_frames)
    return trans_matrices

test_data = np.array([[[1143.52873, 1426.63070],
        [1173.33920, 950.41084],
        [497.55069, 1192.56431],
        [683.19722, 356.49423],
        [962.27561, 953.18146],
        [950.88780, 631.34461],
        [1398.45339, 883.32731],
        [851.62795, 1166.07798],
        [817.24848, 996.53567],
        [655.51762, 1422.98901]],
 [[1184.11036, 971.56517],
        [1154.34418, 1447.65958],
        [508.43994, 1213.77887],
        [693.74653, 377.72513],
        [972.81065, 974.02935],
        [962.01848, 651.73851],
        [1409.44561, 904.52413],
        [862.66761, 1187.06697],
        [828.23598, 1017.57026],
        [621.53039, 939.63127]],
 [[1124.90251, 1443.75824],
        [1154.70004, 967.56144],
        [478.76940, 1209.76671],
        [664.28665, 373.73092],
        [942.82431, 969.96256],
        [932.07174, 648.12686],
        [1379.55540, 900.54206],
        [833.05714, 1183.36170],
        [798.52507, 1013.77275],
        [591.55219, 935.74713]]])

loop_trans_matrices = [get_rigid_transform(i, test_data[0]) for i in test_data]
parellel_trans_matrices = match_frames_joblib_loop(test_data)

np.set_printoptions(suppress=True,
   formatter={'float':'{:.5f}'.format})

print(np.array(loop_trans_matrices[1:]) - np.array(parellel_trans_matrices[1:]))

It gives the output:

[[[-0.00004 0.00004 -0.00628]
  [-0.00004 -0.00004 0.08418]
  [0.00000 0.00000 0.00000]]

 [[-0.00002 -0.00002 0.02500]
  [0.00002 -0.00002 0.01319]
  [0.00000 0.00000 0.00000]]]

I also noticed that the problem disappears if I modify the get_rigid_transform function to include the line aa.NUM_NEAREST_NEIGHBORS = 6.

def get_rigid_transform(frame_sources, reference):
    aa.NUM_NEAREST_NEIGHBORS = 6
    transf, matched_positions = aa.find_transform(frame_sources, reference, max_control_points= 30)
    return transf.params

I don't know if this is an Astroalign issue or something the user should take care of.

dholzber commented 1 year ago

You might want to take a look at #82.