mpound / pdrtpy

PhotoDissociation Region Toolbox Python module
GNU General Public License v3.0
6 stars 3 forks source link

Measurement.read fails if WCS project has 3 dimensions #90

Open mpound opened 10 months ago

mpound commented 10 months ago

reported by Lorenzo Evangelista evangeli@iap.fr

Dear Mark Pound,

I am Pierre Guillard's PhD student working on JWST H2 maps of the CenA galaxy. We’re contacting you via Emilie Habart.

We are trying to generate a pixel-by-pixel gas excitation map by superposing a set of multiple moment maps, however we have an issue with the Measurement.read function.

It looks like the function doesn’t accept our moment maps because they are generated from sliced-cubes with a WCS having 3 dimensions, while 2 are expected. We have tried to modify the Header to make it acceptable but all attempts failed.

Is there a way to allow the function to accept 3 dimension headers?

Thank you very much for your help,

Best wishes,

Lorenzo EVANGELISTA

Here attached you’ll find the map which gives the error. Following you’ll find part of the code and the error.


q = []
nu = []
i=1
for j in intensity:
    infile=f"Moment_maps/Conv_maps/ConvP_S{i}.fits"
    m = Measurement.read(infile,identifier=j).    ######################## Line of error
    q.append(m)
    nu.append(h.upper_colden(m,'cm-2'))
    i+=1
# pick a single scale for all the plots
vmin=q[1].data.min()
vmax=q[1].data.max()

fig,ax = plt.subplots(2,3,subplot_kw={'projection':q[0].wcs},figsize=(9,6))
i = 0
for a in ax:
    for b in a:
        b.imshow(q[i],origin="lower",vmin=vmin,vmax=vmax,cmap='plasma')
        b.set_title(q[i].id)
        i = i+1
plt.rcParams['font.size'] = 8
plt.subplots_adjust(wspace=0.5)

fitted 1 of 1 pixels got 0 exceptions and 0 bad fits WARNING: FITSFixedWarning: The WCS transformation has more axes (3) than the image it is associated with (2) [astropy.wcs.wcs] WARNING: FITSFixedWarning: 'obsfix' made the change 'Set OBSGEO-L to 169.637084 from OBSGEO-[XYZ]. Set OBSGEO-B to 18.461327 from OBSGEO-[XYZ]. Set OBSGEO-H to 1255350594.775 from OBSGEO-[XYZ]'. [astropy.wcs.wcs] Traceback (most recent call last):

File ~/anaconda3/lib/python3.11/site-packages/pdrtpy/measurement.py:650 in fits_measurement_reader z=Measurement(z,unit=z._unit,title=_title)

File ~/anaconda3/lib/python3.11/site-packages/pdrtpy/measurement.py:147 in init self._set_up_for_interp()

File ~/anaconda3/lib/python3.11/site-packages/pdrtpy/measurement.py:338 in _set_up_for_interp self._world_axis = utils.get_xy_from_wcs(self,quantity=False,linear=False)

File ~/anaconda3/lib/python3.11/site-packages/pdrtpy/pdrutils.py:711 in get_xy_from_wcs x=w.array_index_to_world_values(xind,xind)[0]

File ~/anaconda3/lib/python3.11/site-packages/astropy/wcs/wcsapi/low_level_api.py:90 in array_index_to_world_values return self.pixel_to_world_values(*index_arrays[::-1])

File ~/anaconda3/lib/python3.11/site-packages/astropy/wcs/wcsapi/fitswcs.py:336 in pixel_to_world_values world = self.all_pix2world(*pixel_arrays, 0)

File ~/anaconda3/lib/python3.11/site-packages/astropy/wcs/wcs.py:1570 in all_pix2world return self._array_converter(self._all_pix2world, "output", *args, **kwargs)

File ~/anaconda3/lib/python3.11/site-packages/astropy/wcs/wcs.py:1562 in _array_converter raise TypeError(

TypeError: WCS projection has 3 dimensions, so expected 2 (an Nx3 array and the origin argument) or 4 arguments (the position in each dimension, and the origin argument). Instead, 3 arguments were given.

During handling of the above exception, another exception occurred:

Traceback (most recent call last):

File ~/anaconda3/lib/python3.11/site-packages/spyder_kernels/py3compat.py:356 in compat_exec exec(code, globals, locals)

File ~/Desktop/Work/untitled1.py:108 m = Measurement.read(infile,identifier=j)

File ~/anaconda3/lib/python3.11/site-packages/astropy/nddata/mixins/ndio.py:59 in call return self.registry.read(self._cls, *args, **kwargs)

File ~/anaconda3/lib/python3.11/site-packages/astropy/io/registry/core.py:218 in read data = reader(*args, **kwargs)

File ~/anaconda3/lib/python3.11/site-packages/pdrtpy/measurement.py:652 in fits_measurement_reader raise TypeError('could not convert fits_measurement_reader output to Measurement')

TypeError: could not convert fits_measurement_reader output to Measurement

mpound commented 10 months ago

Measurement.read() is failing because your FITS file is not in the appropriate format - it needs 2 HDUs - one for the flux and one for the error. You should use Measurement.make_measurement() first to create the Measurement FITS file from your intensity FITS file. See https://pdrtpy.readthedocs.io/en/latest/pdrtpy.measurement.html#pdrtpy.measurement.Measurement.make_measurement This is also shown in the example notebook: https://github.com/mpound/pdrtpy-nb/blob/master/notebooks/PDRT_Example_Measurements.ipynb in the section "Creating Measurements from FITS files"

Try that and let me know how it goes.

mpound commented 10 months ago

From: guillard@iap.fr Dear Marc

Many thanks for your quick response. I had tried that already, but unfortunately we get the same error: H200S1_flux = utils.get_testdata("H2_0_0_S1_moment0.fits") # H2 0-0 S(1) 17um flux H200S1_combined = "H200S1_flux_error.fits" Measurement.make_measurement(H200S1_flux, error='10%', outfile=H200S1_combined, overwrite=True) H200S1_meas = Measurement.read(H200S1_combined, identifier="H200S1")

Leads to:

File /usr/local/anaconda3/envs/pdrtpy/lib/python3.9/site-packages/astropy/wcs/wcs.py:1365, in WCS.all_pix2world(self, *args, kwargs) 1364 def all_pix2world(self, *args, *kwargs): -> 1365 return self._array_converter( 1366 self._all_pix2world, 'output', args, kwargs)

File /usr/local/anaconda3/envs/pdrtpy/lib/python3.9/site-packages/astropy/wcs/wcs.py:1357, in WCS._array_converter(self, func, sky, ra_dec_order, *args) 1355 return _return_list_of_arrays(axes, origin) -> 1357 raise TypeError( 1358 "WCS projection has {0} dimensions, so expected 2 (an Nx{0} array " 1359 "and the origin argument) or {1} arguments (the position in each " 1360 "dimension, and the origin argument). Instead, {2} arguments were " 1361 "given.".format( 1362 self.naxis, self.naxis + 1, len(args)))

TypeError: WCS projection has 3 dimensions, so expected 2 (an Nx3 array and the origin argument) or 4 arguments (the position in each dimension, and the origin argument). Instead, 3 arguments were given.

Many thanks for your help,

Pierre

mpound commented 10 months ago

example FITS file for testing wcs3.zip