|
| 1 | +import astropy.coordinates as ac |
| 2 | +import astropy.units as au |
| 3 | +import numpy as np |
| 4 | + |
| 5 | + |
| 6 | +def main(): |
| 7 | + lon, lat, height = [], [], [] |
| 8 | + |
| 9 | + with open('brad_verify.txt', 'r') as f: |
| 10 | + for line in f: |
| 11 | + if line.startswith('#'): |
| 12 | + continue |
| 13 | + if line.strip() == '': |
| 14 | + continue |
| 15 | + _, _lon, _lat, _ele = line.split() |
| 16 | + lon.append(float(_lon)) |
| 17 | + lat.append(float(_lat)) |
| 18 | + height.append(float(_ele)) |
| 19 | + |
| 20 | + antennas = ac.EarthLocation.from_geodetic(lon * au.deg, lat * au.deg, height * au.m) |
| 21 | + print(antennas) |
| 22 | + |
| 23 | + solution_file = '/home/albert/Downloads/dsa1650_a_P305_v2.4.6.txt' |
| 24 | + x, y, z = [], [], [] |
| 25 | + with open(solution_file, 'r') as f: |
| 26 | + for line in f: |
| 27 | + if line.startswith('#'): |
| 28 | + continue |
| 29 | + if line.strip() == '': |
| 30 | + continue |
| 31 | + _x, _y, _z = line.split(',') |
| 32 | + x.append(float(_x)) |
| 33 | + y.append(float(_y)) |
| 34 | + z.append(float(_z)) |
| 35 | + |
| 36 | + antennas_sol = ac.EarthLocation.from_geocentric(x * au.m, y * au.m, z * au.m) |
| 37 | + print(antennas_sol) |
| 38 | + |
| 39 | + # For each antenna in `antennas`, find the closest antenna in `antennas_sol` |
| 40 | + antennas_itrs = antennas.get_itrs().cartesian.xyz.T |
| 41 | + antennas_sol_itrs = antennas_sol.get_itrs().cartesian.xyz.T |
| 42 | + dist = [] |
| 43 | + for i, antenna in enumerate(antennas_itrs): |
| 44 | + # Calculate the distance to each antenna in antennas_sol |
| 45 | + distances = np.linalg.norm(antennas_sol_itrs - antenna, axis=1) |
| 46 | + # Find the index of the closest antenna |
| 47 | + closest_index = np.argmin(distances) |
| 48 | + # Print the result |
| 49 | + print( |
| 50 | + f"Antenna {i}: Closest antenna in solution is {closest_index} with distance {distances[closest_index]:.2f}") |
| 51 | + dist.append(distances[closest_index].value) |
| 52 | + |
| 53 | + # Print the average distance |
| 54 | + print(f"Average distance: {np.mean(dist):.2f}") |
| 55 | + # Print the maximum distance |
| 56 | + print(f"Maximum distance: {np.max(dist):.2f}") |
| 57 | + # Print the minimum distance |
| 58 | + print(f"Minimum distance: {np.min(dist):.2f}") |
| 59 | + |
| 60 | + |
| 61 | +if __name__ == '__main__': |
| 62 | + main() |
0 commit comments