Now using format digits from header template (45e4f63e) · Commits · noctis / noche

noche/noche.py

+91 −56

Original line number	Diff line number	Diff line
		@@ -50,6 +50,9 @@ class Noche:

		if debug:
		log.setLevel('DEBUG')
		else:
		log.setLevel('INFO')


		log.debug(sys._getframe().f_code.co_name)

		@@ -60,6 +63,8 @@ class Noche:
		self._head_dir = "headers"
		self._obs_dir = "observatories"

		self._header_dict = {}

		self.header = fits.Header()
		self.hdu = None

		@@ -95,7 +100,7 @@ class Noche:
		log.debug(sys._getframe().f_code.co_name)

		if not path:
		path = Path(__file__).parent / self._head_dir / "header_base_v1.ini"
		path = Path(__file__).parent / self._head_dir / "header_base_v1_with_format.ini"
		log.info(path)

		config = self._load_config(path)
		@@ -103,8 +108,8 @@ class Noche:
		for section in config.sections():
		for k, value in config.items(section):
		try:
		val, comment = value.split("\|")
		# val, typ, fmt, comment = value.split("\|")
		# val, comment = value.split("\|")
		val, typ, fmt, comment = value.split("\|")
		comment = comment.strip()
		except ValueError:
		# HISTORY and COMMENT have no comment
		@@ -113,10 +118,43 @@ class Noche:
		val = self._parse(val)
		log.debug(f"{k:<11}: initializing with: {val} / {comment}")
		self.header[k] = val, comment
		# Filling the header dict with format info too.
		self._header_dict[k] = (val, typ, fmt, comment)

		self._update()


		def fill_keyword(self, k, val):
		"""
		Fill new header keyword with its value.
		Why use this function instea of simply self.header[keyword] = value?
		Because it also takes care of the format, rounding,
		and and significant digits
		taken from header tempalte.

		Parameters
		----------
		k : str
		Header keyword
		val : str or bool or int or float
		value to insert in the fits header
		"""

		# Getting from the dict the type and the digits.
		# default val (_) and comment( __) are not used.
		_, typ, fmt, __ = self._header_dict[k.lower()]

		val = self._parse(val)

		# log.debug(f"{val}, {typ}, {fmt}")

		if typ.strip() == "float":
		val = header_round(val, int(fmt))

		self.header[k] = val
		log.info(f"{k:<11}: filled with {val} as {typ:<5}, {fmt:<1} digits")


		def load_noctis_observatory(self, name='oarpaf', fits_file=None):
		"""
		Load one of the NOCTIS observatory parameters such as
		@@ -170,14 +208,14 @@ class Noche:
		log.error(f"{k} not in header template")
		continue

		val = self._parse(loc[k])
		val = loc[k]

		if val is None:
		log.warning(f"{k:<8} : fixed value missing! Leave blank")
		pass
		else:
		log.info(f"{k:<11} : set to fixed value, value: {val:<28}")
		self.header[k] = val
		self.fill_keyword(k, val)

		if fits_file:
		self.fill_from_fits_file(path, fits_file)
		@@ -223,9 +261,8 @@ class Noche:
		fits_keyword = loc[k]

		try:
		fits_value = fits_file_header[fits_keyword]
		val = self._parse(fits_value)
		self.header[k] = val
		val = fits_file_header[fits_keyword]
		self.fill_keyword(k, val)
		log.info(f"{k:<11} : mapped to {fits_keyword:<11}, value: {val:<28}")
		except TypeError as e:
		log.debug(e)
		@@ -239,11 +276,10 @@ class Noche:
		pass

		if k in pre_formulas:
		x = self.header[k]
		digits = len(str(x).split(".")[1])
		x = self.header[k] # i.e. value of RA: x=12.345
		new_val = eval(pre_formulas[k])
		try:
		self.header[k] = header_round(new_val, digits)
		self.fill_keyword(k, new_val)
		log.warning(f"{k:<11}: pre-formula from {x:<28} to {self.header[k]:<28}")
		except Exception as e:
		log.error(f"Failed to evaluate pre-formula for {k}: {e}")
		@@ -256,17 +292,16 @@ class Noche:
		for k, formula in post_formulas.items():
		if k in self.header:
		x = self.header[k]
		digits = len(str(x).split(".")[1])
		new_val = eval(formula)
		try:
		self.header[k] = header_round(new_val, digits)
		self.fill_keyword(k, new_val)
		log.warning(f"{k:<8}: tweak from {x:<28} to {self.header[k]:<28}")
		except Exception as e:
		log.error(f"Failed to evaluate tweak for {k}: {e}")
		else:
		log.warning(f"{k:<8}: skipped as not in header")

		self.header["FILEORIG"] = filename.name
		self.fill_keyword("FILEORIG", filename.name)


		def set_location(self, lon, lat, alt):
		@@ -287,9 +322,9 @@ class Noche:

		self._location = EarthLocation(lon, lat, alt)

		self.header["OBS-LONG"] = lon
		self.header["OBS-LAT"] = lat
		self.header["OBS-ELEV"] = alt
		self.fill_keyword("OBS-LONG", lon)
		self.fill_keyword("OBS-LAT", lat)
		self.fill_keyword("OBS-ELEV", alt)


		def set_obstime(self, obstime):
		@@ -308,9 +343,9 @@ class Noche:
		if self._coord != None:
		self._coord.obstime = time

		self.header['DATE'] = time.isot.split("T")[0]
		self.header['DATE-OBS'] = time.isot
		self.header['MJD-OBS'] = time.mjd
		self.fill_keyword("DATE", time.isot.split("T")[0])
		self.fill_keyword("DATE-OBS", time.isot)
		self.fill_keyword("MJD-OBS", time.mjd)

		#self._update()

		@@ -352,7 +387,7 @@ class Noche:
		log.error("Cannot resolve name")
		return [None, None]

		self.header["OBJECT"] = objname
		self.fill_keyword("OBJECT", objname)


		def set_coordinates(self, ra, dec, obstime=None):
		@@ -380,12 +415,12 @@ class Noche:

		with custom_float():

		self.header['RA'] = coord.ra.to_string(unit=u.hourangle, sep=':',
		pad=True, precision=1)
		self.header['DEC'] = coord.dec.to_string(unit=u.deg, sep=':',
		pad=True, precision=1)
		self.header['RA_DEG'] = header_round(coord.ra.deg, 7)
		self.header['DEC_DEG'] = header_round(coord.dec.deg,7)
		self.fill_keyword("RA", coord.ra.to_string(unit=u.hourangle, sep=':',
		pad=True, precision=1) )
		self.fill_keyword("DEC", coord.dec.to_string(unit=u.deg, sep=':',
		pad=True, precision=1) )
		self.fill_keyword("RA_DEG", coord.ra.deg)
		self.fill_keyword("DEC_DEG", coord.dec.deg)

		self._update()

		@@ -411,29 +446,29 @@ class Noche:
		# Altitudine and Azimuth
		with custom_float():

		self.header['ALT'] = header_round(altaz.alt.deg, 7)
		self.header['AZ'] = header_round(altaz.az.deg, 7)
		self.header['AIRMASS'] = header_round(altaz.secz.value, 2)
		self.fill_keyword("ALT", altaz.alt.deg)
		self.fill_keyword("AZ", altaz.az.deg)
		self.fill_keyword("AIRMASS", altaz.secz.value)

		# Local Sideral Time and Hour Angle
		lst = self._obstime.sidereal_time('mean', longitude=self._location.lon)
		ha = (lst - self._coord.ra).hour
		lst_hours = lst.hour

		self.header['LST'] = header_round(lst_hours, 2)
		self.header['HA'] = header_round(ha, 2)
		self.fill_keyword("LST", lst_hours)
		self.fill_keyword("HA", ha)

		# Position angle: with respect to Celestial North Pole
		north_celestial = SkyCoord(ra=0u.deg, dec=90u.deg, frame='icrs')
		posangle = self._coord.position_angle(north_celestial).to(u.deg).value
		self.header['POSANGLE'] = header_round(posangle, 2)
		self.fill_keyword("POSANGLE", posangle)

		# Parallactic angle: between local meridian and celestial axis
		parangle = (posangle - altaz.az.deg + 360) % 360
		if parangle > 180:
		parangle -= 360 # Wrap to [-180, 180]

		self.header['PARANGLE'] = header_round(parangle, 2)
		self.fill_keyword("PARANGLE", parangle)


		def set_wcs(self, angle=None):
		@@ -455,17 +490,17 @@ class Noche:
		if self._coord == None or self._location == None:
		raise ValueError("Observation Coordinates, Instrument parameters must be set.")

		detsize = self.header['DETSIZE'].strip('[]')
		x_str, y_str = detsize.split(',')
		detsize = self.header["DETSIZE"].strip("[]")
		x_str, y_str = detsize.split(",")

		# Ottieni i limiti numerici
		_, xsize = map(int, x_str.split(':'))
		_, ysize = map(int, y_str.split(':'))
		_, xsize = map(int, x_str.split(":"))
		_, ysize = map(int, y_str.split(":"))

		with custom_float():
		crpix = [xsize/self.header["XBINNING"]/2, ysize/self.header["YBINNING"]/2]
		cdelt1 = header_round(self.header["PIXSCALE"]self.header["XBINNING"]u.arcsec.to(u.deg), 7)
		cdelt2 = header_round(self.header["PIXSCALE"]self.header["YBINNING"]u.arcsec.to(u.deg), 7)
		cdelt1 = self.header["PIXSCALE"]self.header["XBINNING"]u.arcsec.to(u.deg)
		cdelt2 = self.header["PIXSCALE"]self.header["YBINNING"]u.arcsec.to(u.deg)

		if not angle:
		angle = self.header["DEROTANG"] + self.header["DETROT"]
		@@ -481,16 +516,16 @@ class Noche:
		flip = -1 # East to the left

		with custom_float():
		self.header['CRPIX1'] = header_round(crpix[0], 7)
		self.header['CRPIX2'] = header_round(crpix[1], 7)
		self.header['CRVAL1'] = header_round(crval_ra, 7)
		self.header['CRVAL2'] = header_round(crval_dec, 7)
		self.header['CDELT1'] = header_round(cdelt1 * flip, 7)
		self.header['CDELT2'] = header_round(cdelt2, 7)
		self.header["PC1_1"] = header_round(+np.cos(angle), 5)
		self.header["PC1_2"] = header_round(-np.sin(angle), 5)
		self.header["PC2_1"] = header_round(+np.sin(angle), 5)
		self.header["PC2_2"] = header_round(+np.cos(angle), 5)
		self.fill_keyword("CRPIX1", crpix[0])
		self.fill_keyword("CRPIX2", crpix[1])
		self.fill_keyword("CRVAL1", crval_ra)
		self.fill_keyword("CRVAL2", crval_dec)
		self.fill_keyword("CDELT1", cdelt1 * flip)
		self.fill_keyword("CDELT2", cdelt2)
		self.fill_keyword("PC1_1", +np.cos(angle))
		self.fill_keyword("PC1_2", -np.sin(angle))
		self.fill_keyword("PC2_1", +np.sin(angle))
		self.fill_keyword("PC2_2", +np.cos(angle))


		def set_ambient(self):
		@@ -522,16 +557,16 @@ class Noche:

		# MOONDIST: angular distance between target and Moon
		moondist = moon.separation(self._coord, origin_mismatch="ignore").deg
		self.header['MOONDIST'] = header_round(moondist, 1)
		self.fill_keyword("MOONDIST", moondist)

		# MOONPHAS: Moon phase
		elongation = moon.separation(sun).deg
		moonphas = (1 + np.cos(np.radians(elongation))) / 2
		self.header['MOONPHAS'] = header_round(moonphas, 2)
		self.fill_keyword("MOONPHAS", moonphas)

		# SUNALT: Sun altitude above the horizon
		sunalt = sun.altaz.alt.deg
		self.header['SUNALT'] = header_round(sunalt, 1)
		self.fill_keyword("SUNALT", sunalt)


		def check_empty(self):