Skip to content

Commit

Permalink
repair make starfield
Browse files Browse the repository at this point in the history
  • Loading branch information
Jérémy Neveu committed Jan 28, 2024
1 parent ea07482 commit e9bf965
Show file tree
Hide file tree
Showing 2 changed files with 222 additions and 52 deletions.
94 changes: 42 additions & 52 deletions spectractor/simulation/image_simulation.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -122,64 +122,54 @@ def __init__(self, base_image, flux_factor=1):
def set_star_list(self):
x0, y0 = self.image.target_pixcoords
sources_file_name = set_sources_file_name(self.image.file_name)
if os.path.isfile(sources_file_name):
wcs_file_name = set_wcs_file_name(self.image.file_name)
gaia_catalog_file_name = set_gaia_catalog_file_name(self.image.file_name)
if os.path.isfile(wcs_file_name) and os.path.isfile(gaia_catalog_file_name):
# load gaia catalog
gaia_catalog = ascii.read(gaia_catalog_file_name, format="ecsv")
gaia_coord_after_motion = get_gaia_coords_after_proper_motion(gaia_catalog, self.image.date_obs)
# load WCS
wcs = load_wcs_from_file(wcs_file_name)
# catalog matching to set star positions using Gaia
target_coord = wcs.all_pix2world([x0], [y0], 0)
target_coord = SkyCoord(ra=target_coord[0] * units.deg, dec=target_coord[1] * units.deg,
frame="icrs", obstime=self.image.date_obs, equinox="J2000")
gaia_target_index, dist_2d, dist_3d = target_coord.match_to_catalog_sky(gaia_coord_after_motion)
dx, dy = 0, 0
for gaia_i in range(len(gaia_catalog)):
x, y = wcs.all_world2pix(gaia_coord_after_motion[gaia_i].ra,
gaia_coord_after_motion[gaia_i].dec, 0)
if gaia_i == gaia_target_index[0]:
dx = x0 - x
dy = y0 - y
A = 10 ** (-gaia_catalog['phot_g_mean_mag'][gaia_i] / 2.5)
self.stars.append(StarModel([x, y], self.image.target_star2D, A))
self.pixcoords.append([x, y])
# rescale using target fitted amplitude
amplitudes = np.array([star.amplitude for star in self.stars])
target_flux = self.image.target_star2D.params.values[0]
amplitudes *= target_flux / self.stars[gaia_target_index[0]].amplitude * self.flux_factor
for k, star in enumerate(self.stars):
star.amplitude = amplitudes[k]
# shift x,y star positions according to target position
star.x0 += dx
star.y0 += dy
star.psf.params.values[1] += dx
star.psf.params.values[2] += dy
star.psf.params.values[0] = amplitudes[k]
elif os.path.isfile(sources_file_name):
# load sources positions and flux
sources = Table.read(sources_file_name)
sources['X'].name = "xcentroid"
sources['Y'].name = "ycentroid"
sources['FLUX'].name = "flux"
# test presence of WCS and gaia catalog files
wcs_file_name = set_wcs_file_name(self.image.file_name)
gaia_catalog_file_name = set_gaia_catalog_file_name(self.image.file_name)
if os.path.isfile(wcs_file_name) and os.path.isfile(gaia_catalog_file_name):
# load gaia catalog
gaia_catalog = ascii.read(gaia_catalog_file_name, format="ecsv")
gaia_coord_after_motion = get_gaia_coords_after_proper_motion(gaia_catalog, self.image.date_obs)
# load WCS
wcs = load_wcs_from_file(wcs_file_name)
# catalog matching to set star positions using Gaia
sources_coord = wcs.all_pix2world(sources['xcentroid'], sources['ycentroid'], 0)
target_coord = wcs.all_pix2world([x0], [y0], 0)
sources_coord = SkyCoord(ra=sources_coord[0] * units.deg, dec=sources_coord[1] * units.deg,
frame="icrs", obstime=self.image.date_obs, equinox="J2000")
target_coord = SkyCoord(ra=target_coord[0] * units.deg, dec=target_coord[1] * units.deg,
frame="icrs", obstime=self.image.date_obs, equinox="J2000")
gaia_index, dist_2d, dist_3d = sources_coord.match_to_catalog_sky(gaia_coord_after_motion)
gaia_target_index, dist_2d, dist_3d = target_coord.match_to_catalog_sky(gaia_coord_after_motion)
#for k, gaia_i in enumerate(gaia_index):
# x, y = wcs.all_world2pix(gaia_coord_after_motion[gaia_i].ra, gaia_coord_after_motion[gaia_i].dec, 0)
# A = sources['flux'][k] * self.flux_factor
# self.stars.append(StarModel([x, y], self.image.target_star2D, A))
# self.pixcoords.append([x, y])
dx, dy = 0, 0
for gaia_i in range(len(gaia_catalog)):
x, y = wcs.all_world2pix(gaia_coord_after_motion[gaia_i].ra,
gaia_coord_after_motion[gaia_i].dec, 0)
if gaia_i == gaia_target_index[0]:
dx = x0 - x
dy = y0 - y
A = 10**(-gaia_catalog['phot_g_mean_mag'][gaia_i]/2.5)
self.stars.append(StarModel([x, y], self.image.target_star2D, A))
self.pixcoords.append([x, y])
# rescale using target fitted amplitude
amplitudes = np.array([star.amplitude for star in self.stars])
target_flux = self.image.target_star2D.params.values[0]
amplitudes *= target_flux / self.stars[gaia_target_index[0]].amplitude * self.flux_factor
for k, star in enumerate(self.stars):
star.amplitude = amplitudes[k]
star.x0 += dx
star.y0 += dy
star.psf.params.values[1] += dx
star.psf.params.values[2] += dy
star.psf.params.values[0] = amplitudes[k]

else:
for k, source in enumerate(sources):
x, y = sources['xcentroid'][k], sources['ycentroid'][k]
A = sources['flux'][k] * self.flux_factor
self.stars.append(StarModel([x, y], self.image.target_star2D, A))
self.pixcoords.append([x, y])
for k, source in enumerate(sources):
x, y = sources['xcentroid'][k], sources['ycentroid'][k]
A = sources['flux'][k] * self.flux_factor
self.stars.append(StarModel([x, y], self.image.target_star2D, A))
self.pixcoords.append([x, y])
else:
# try extraction using iraf source detection
# mask background, faint stars, and saturated pixels
data = np.copy(self.image.data)
# mask order0 and spectrum
Expand Down
Loading

0 comments on commit e9bf965

Please sign in to comment.