Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

The print(wcs.tofits_sip()) failed in your example in https://gwcs.readthedocs.io/en/latest/#a-step-by-step-example-of-constructing-an-imaging-gwcs-object #445

Open
lznakano opened this issue Mar 13, 2023 · 2 comments

Comments

@lznakano
Copy link

I followed the example in https://gwcs.readthedocs.io/en/latest/#a-step-by-step-example-of-constructing-an-imaging-gwcs-object to make a wcs. When I tried to print it out, it failed with a lot of errors.
Here is the function I defined (copied and pasted from the example):

from astropy.io import fits
import astropy.time as time
import os
from astropy.modeling import models
import asdf
from astropy import units as u
import astropy.time as atime
import roman_datamodels.datamodels as rdm
import roman_datamodels.stnode as stnode
import numpy as np
from astropy import coordinates as coord
from gwcs import wcs
from gwcs import coordinate_frames as cf
from gwcs.utils import make_fitswcs_transform
from roman_datamodels import units as ru

#compose a wcs object
def make_wcs():
   shift_by_crpix = models.Shift(-(2048 - 1)*u.pix) & models.Shift(-(1024 - 1)*u.pix)
   matrix = np.array([[1.290551569736E-05, 5.9525007864732E-06],
                   [5.0226382102765E-06 , -1.2644844123757E-05]])
   rotation = models.AffineTransformation2D(matrix * u.deg,
                                         translation=[0, 0] * u.deg)
   rotation.input_units_equivalencies = {"x": u.pixel_scale(1*u.deg/u.pix),
                                      "y": u.pixel_scale(1*u.deg/u.pix)}
   rotation.inverse = models.AffineTransformation2D(np.linalg.inv(matrix) * u.pix,
                                                 translation=[0, 0] * u.pix)
   rotation.inverse.input_units_equivalencies = {"x": u.pixel_scale(1*u.pix/u.deg),
                                              "y": u.pixel_scale(1*u.pix/u.deg)}

   tan = models.Pix2Sky_TAN()
   celestial_rotation =  models.RotateNative2Celestial(5.63056810618*u.deg, -72.05457184279*u.deg, 180*u.deg)
   det2sky = shift_by_crpix | rotation | tan | celestial_rotation
   det2sky.name = "linear_transform"


   detector_frame = cf.Frame2D(name="detector", axes_names=("x", "y"),
                            unit=(u.pix, u.pix))
   sky_frame = cf.CelestialFrame(reference_frame=coord.ICRS(), name='icrs',
                              unit=(u.deg, u.deg))

   pipeline = [(detector_frame, det2sky),
            (sky_frame, None)
           ]
   wcsobj = wcs.WCS(pipeline)
   wcsobj.bounding_box = ((0, 2048), (0,1024))
   return wcsobj`

wcs=make_wcs()
print(wcs.to_fits_sip())

File ~/miniconda3/envs/roman-data-workshop-env/lib/python3.10/site-packages/gwcs/wcs.py:1609, in WCS.to_fits_sip(self, bounding_box, max_pix_error, degree, max_inv_pix_error, inv_degree, npoints, crpix, projection, verbose)
1606 if celestial_group is None:
1607 raise ValueError("The to_fits_sip requires an output celestial frame.")
-> 1609 hdr = self._to_fits_sip(
1610 celestial_group=celestial_group,
1611 keep_axis_position=False,
1612 bounding_box=bounding_box,
1613 max_pix_error=max_pix_error,
1614 degree=degree,
1615 max_inv_pix_error=max_inv_pix_error,
1616 inv_degree=inv_degree,
1617 npoints=npoints,
1618 crpix=crpix,
1619 projection=projection,
1620 matrix_type='CD',
1621 verbose=verbose
1622 )
1624 return hdr

File ~/miniconda3/envs/roman-data-workshop-env/lib/python3.10/site-packages/gwcs/wcs.py:1797, in WCS._to_fits_sip(self, celestial_group, keep_axis_position, bounding_box, max_pix_error, degree, max_inv_pix_error, inv_degree, npoints, crpix, projection, matrix_type, verbose)
1794 if (xmax - xmin) < 1 or (ymax - ymin) < 1:
1795 raise ValueError("Bounding box is too small for fitting a SIP polynomial")
-> 1797 lon, lat = transform(crpix1, crpix2)
1799 # Now rotate to native system and deproject. Recall that transform
1800 # expects pixels in the original coordinate system, but the SIP
1801 # transform is relative to crpix coordinates, thus the initial shift.
1802 ntransform = ((Shift(crpix1) & Shift(crpix2)) | transform
1803 | RotateCelestial2Native(lon, lat, 180)
1804 | sky2pix_proj)

File ~/miniconda3/envs/roman-data-workshop-env/lib/python3.10/site-packages/astropy/modeling/core.py:1136, in Model.call(self, *args, **kwargs)
1131 # prepare for model evaluation (overridden in CompoundModel)
1132 evaluate, inputs, broadcasted_shapes, kwargs = self._pre_evaluate(
1133 *args, **kwargs
1134 )
-> 1136 outputs = self._generic_evaluate(evaluate, inputs, fill_value, with_bbox)
1138 # post-process evaluation results (overridden in CompoundModel)
1139 return self._post_evaluate(
1140 inputs, outputs, broadcasted_shapes, with_bbox, **kwargs
1141 )

File ~/miniconda3/envs/roman-data-workshop-env/lib/python3.10/site-packages/astropy/modeling/core.py:1098, in Model._generic_evaluate(self, evaluate, _inputs, fill_value, with_bbox)
1096 outputs = bbox.evaluate(evaluate, _inputs, fill_value)
1097 else:
-> 1098 outputs = evaluate(_inputs)
1099 return outputs

File ~/miniconda3/envs/roman-data-workshop-env/lib/python3.10/site-packages/astropy/modeling/core.py:3352, in CompoundModel._pre_evaluate..evaluate(_inputs)
3351 def evaluate(_inputs):
-> 3352 return self._evaluate(*_inputs, **kwargs)

File ~/miniconda3/envs/roman-data-workshop-env/lib/python3.10/site-packages/astropy/modeling/core.py:3380, in CompoundModel._evaluate(self, *args, **kw)
3378 if op != "fix_inputs":
3379 if op != "&":
-> 3380 leftval = self.left(*args, **kw)
3381 if op != "|":
3382 rightval = self.right(*args, **kw)

File ~/miniconda3/envs/roman-data-workshop-env/lib/python3.10/site-packages/astropy/modeling/core.py:1136, in Model.call(self, *args, **kwargs)
1131 # prepare for model evaluation (overridden in CompoundModel)
1132 evaluate, inputs, broadcasted_shapes, kwargs = self._pre_evaluate(
1133 *args, **kwargs
1134 )
-> 1136 outputs = self._generic_evaluate(evaluate, inputs, fill_value, with_bbox)
1138 # post-process evaluation results (overridden in CompoundModel)
1139 return self._post_evaluate(
1140 inputs, outputs, broadcasted_shapes, with_bbox, **kwargs
1141 )

File ~/miniconda3/envs/roman-data-workshop-env/lib/python3.10/site-packages/astropy/modeling/core.py:1098, in Model._generic_evaluate(self, evaluate, _inputs, fill_value, with_bbox)
1096 outputs = bbox.evaluate(evaluate, _inputs, fill_value)
1097 else:
-> 1098 outputs = evaluate(_inputs)
1099 return outputs

File ~/miniconda3/envs/roman-data-workshop-env/lib/python3.10/site-packages/astropy/modeling/core.py:3352, in CompoundModel._pre_evaluate..evaluate(_inputs)
3351 def evaluate(_inputs):
-> 3352 return self._evaluate(*_inputs, **kwargs)

File ~/miniconda3/envs/roman-data-workshop-env/lib/python3.10/site-packages/astropy/modeling/core.py:3380, in CompoundModel._evaluate(self, *args, **kw)
3378 if op != "fix_inputs":
3379 if op != "&":
-> 3380 leftval = self.left(*args, **kw)
3381 if op != "|":
3382 rightval = self.right(*args, **kw)

[... skipping similar frames: Model.__call__ at line 1136 (2 times), Model._generic_evaluate at line 1098 (2 times), CompoundModel._pre_evaluate.<locals>.evaluate at line 3352 (2 times), CompoundModel._evaluate at line 3380 (1 times)]

File ~/miniconda3/envs/roman-data-workshop-env/lib/python3.10/site-packages/astropy/modeling/core.py:3380, in CompoundModel._evaluate(self, *args, **kw)
3378 if op != "fix_inputs":
3379 if op != "&":
-> 3380 leftval = self.left(*args, **kw)
3381 if op != "|":
3382 rightval = self.right(*args, **kw)

File ~/miniconda3/envs/roman-data-workshop-env/lib/python3.10/site-packages/astropy/modeling/core.py:1136, in Model.call(self, *args, **kwargs)
1131 # prepare for model evaluation (overridden in CompoundModel)
1132 evaluate, inputs, broadcasted_shapes, kwargs = self._pre_evaluate(
1133 *args, **kwargs
1134 )
-> 1136 outputs = self._generic_evaluate(evaluate, inputs, fill_value, with_bbox)
1138 # post-process evaluation results (overridden in CompoundModel)
1139 return self._post_evaluate(
1140 inputs, outputs, broadcasted_shapes, with_bbox, **kwargs
1141 )

File ~/miniconda3/envs/roman-data-workshop-env/lib/python3.10/site-packages/astropy/modeling/core.py:1098, in Model._generic_evaluate(self, evaluate, _inputs, fill_value, with_bbox)
1096 outputs = bbox.evaluate(evaluate, _inputs, fill_value)
1097 else:
-> 1098 outputs = evaluate(_inputs)
1099 return outputs

File ~/miniconda3/envs/roman-data-workshop-env/lib/python3.10/site-packages/astropy/modeling/core.py:3352, in CompoundModel._pre_evaluate..evaluate(_inputs)
3351 def evaluate(_inputs):
-> 3352 return self._evaluate(*_inputs, **kwargs)

File ~/miniconda3/envs/roman-data-workshop-env/lib/python3.10/site-packages/astropy/modeling/core.py:3387, in CompoundModel._evaluate(self, args, **kw)
3384 rightval = None
3386 else:
-> 3387 leftval = self.left(
(args[: self.left.n_inputs]), **kw)
3388 rightval = self.right(*(args[self.left.n_inputs :]), **kw)
3390 if op != "|":

File ~/miniconda3/envs/roman-data-workshop-env/lib/python3.10/site-packages/astropy/modeling/core.py:424, in call(self, model_set_axis, with_bounding_box, fill_value, equivalencies, inputs_map, *inputs, **new_inputs)
415 args = ("self",)
416 kwargs = {
417 "model_set_axis": None,
418 "with_bounding_box": False,
(...)
421 "inputs_map": None,
422 }
--> 424 new_call = make_function_with_signature(
425 call, args, kwargs, varargs="inputs", varkwargs="new_inputs"
426 )
428 # The following makes it look like call
429 # was defined in the class
430 update_wrapper(new_call, cls)

File ~/miniconda3/envs/roman-data-workshop-env/lib/python3.10/site-packages/astropy/modeling/core.py:401, in _ModelMeta._handle_special_methods..call(self, *inputs, **kwargs)
399 def call(self, *inputs, **kwargs):
400 """Evaluate this model on the supplied inputs."""
--> 401 return super(cls, self).call(*inputs, **kwargs)

File ~/miniconda3/envs/roman-data-workshop-env/lib/python3.10/site-packages/astropy/modeling/core.py:1132, in Model.call(self, *args, **kwargs)
1129 fill_value = kwargs.pop("fill_value", np.nan)
1131 # prepare for model evaluation (overridden in CompoundModel)
-> 1132 evaluate, inputs, broadcasted_shapes, kwargs = self._pre_evaluate(
1133 *args, **kwargs
1134 )
1136 outputs = self._generic_evaluate(evaluate, inputs, fill_value, with_bbox)
1138 # post-process evaluation results (overridden in CompoundModel)

File ~/miniconda3/envs/roman-data-workshop-env/lib/python3.10/site-packages/astropy/modeling/core.py:978, in Model._pre_evaluate(self, *args, **kwargs)
973 """
974 Model specific input setup that needs to occur prior to model evaluation
975 """
977 # Broadcast inputs into common size
--> 978 inputs, broadcasted_shapes = self.prepare_inputs(*args, **kwargs)
980 # Setup actual model evaluation method
981 parameters = self._param_sets(raw=True, units=True)

File ~/miniconda3/envs/roman-data-workshop-env/lib/python3.10/site-packages/astropy/modeling/core.py:2126, in Model.prepare_inputs(self, model_set_axis, equivalencies, *inputs, **kwargs)
2122 self._validate_input_shapes(inputs, self.inputs, model_set_axis)
2124 inputs_map = kwargs.get("inputs_map", None)
-> 2126 inputs = self._validate_input_units(inputs, equivalencies, inputs_map)
2128 # The input formatting required for single models versus a multiple
2129 # model set are different enough that they've been split into separate
2130 # subroutines
2131 if self._n_models == 1:

File ~/miniconda3/envs/roman-data-workshop-env/lib/python3.10/site-packages/astropy/modeling/core.py:2222, in Model._validate_input_units(self, inputs, equivalencies, inputs_map)
2216 if (
2217 not self.input_units_allow_dimensionless[input_name]
2218 and input_unit is not dimensionless_unscaled
2219 and input_unit is not None
2220 ):
2221 if np.any(inputs[i] != 0):
-> 2222 raise UnitsError(
2223 f"{name}: Units of input '{self.inputs[i]}',"
2224 " (dimensionless), could not be converted to required "
2225 f"input units of {input_unit} "
2226 f"({input_unit.physical_type})"
2227 )
2228 return inputs

UnitsError: Shift: Units of input 'x', (dimensionless), could not be converted to required input units of pix (unknown)

@nden
Copy link
Collaborator

nden commented Mar 21, 2023

@lznakano This is a bug (or feature) of astropy.models and we are looking into work arounds or fixes.

@lznakano
Copy link
Author

lznakano commented Mar 21, 2023 via email

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

2 participants