Open kecnry opened 1 year ago
I don't understand why the operator is so not caring about the type of other
that is allowed for proper math. self.add()
here is all the way in astropy.nddata
and it has no knowledge of what a spectral axis is, so any spectral axis check has to happen before we call self.add()
. But it is not clear to me why you accept seemingly any random kind of data. How am I supposed to know where the spectral axis is in any given kind of things accepted here? The arithmetric tests only testing Spectrum1D
inputs.
Also, this is out of scope, right?
My proposal for __add__
and __sub__
:
other
is not Spectrum1D
.sys.float_info.epsilon
good here?nddata
but I'll have to double check.Is __mul__
and __div__
in scope here? The original post above only mentioned -
operator. Those will have very different logic from __add__
and __sub__
.
sr
is okay if it doesn't mess with the rest of the specutils
calculations). Unitless scalar should be acceptable. Otherwise, it must be Spectrum1D
and spectral axis must match as above.self
. The result here would give you a unitless Spectrum1D
. Do you support unitless Spectrum1D
at all? I am not talking about count
. It has to be u.dimensionless_unscaled
.Do I have to worry about magnitude flux units in specutils
?
Thanks for writing up a proposal @pllim. Here are my thoughts:
My proposal for
__add__
and__sub__
:* Throw error if `other` is not `Spectrum1D`.
I think we should also allow other
to be a Quantity
of matching unit/shape to the flux
array. I'm thinking about the case where someone wants to do continuum-subtraction - I don't think the user should necessarily have to put their continuum into a Spectrum1D
after average a flux array or something. If the other
operand is also a Spectrum1D
, then we should check that the axes match.
* As Ricky suggested, strict check for spectral axis values (match must be 1-to-1 within tolerance). Is `sys.float_info.epsilon` good here?
I haven't used this before, from a quick look it's probably sufficient.
Is
__mul__
and__div__
in scope here? The original post above only mentioned-
operator. Those will have very different logic from__add__
and__sub__
.
These are in scope, but the unit compatibility checks you're proposing below aren't. If you want to add those now, go for it, but the scope of this issue is checking that the spectral axes match.
* Other must be unitless or equivalent (e.g., maybe `sr` is okay if it doesn't mess with the rest of the `specutils` calculations). Unitless scalar should be acceptable. Otherwise, it must be `Spectrum1D` and spectral axis must match as above. * For division only, if other is not unitless, then it must have the same unit as `self`. The result here would give you a unitless `Spectrum1D`. Do you support unitless `Spectrum1D` at all? I am not talking about `count`. It has to be `u.dimensionless_unscaled`.
We do use u.dimensionless_unscaled
in a few places.
Do I have to worry about magnitude flux units in
specutils
?
No.
I think we should also allow other to be a Quantity of matching unit/shape
Can you please clarify? It makes sense to allow scalar Quantity, but if the Quantity has an array, it feels wrong to just assume user knows the arrays maps to the same spectral axis.
the unit compatibility checks you're proposing below aren't
Well, I guess they can also be addressed in a different PR...
I think we should also allow other to be a Quantity of matching unit/shape
Can you please clarify? It makes sense to allow scalar Quantity, but if the Quantity has an array, it feels wrong to just assume user knows the arrays maps to the same spectral axis.
Here's my theoretical workflow: I have a 2D or 3D Spectrum1D
, I take an average of some background region, I create a numpy array of the same shape as flux
and assign that average to every array element and turn it into a Quantity
using the flux
unit. Now I want to subtract it from the original spectrum - I know it's something that I want to apply to every element of the spectrum and might be annoyed that I have to take the extra step of getting the spectral axis from the original to turn this into a Spectrum1D
just to subtract it.
It seems unlikely (thought not impossible) that I would just coincidentally have something of the exact same shape as my original spectrum and be bamboozled by addition/subtraction working when it shouldn't.
the unit compatibility checks you're proposing below aren't
Well, I guess they can also be addressed in a different PR...
Right, I'm not saying it shouldn't be done, I'm just saying if you want to minimize the scope of your PR for this issue, the units stuff can be separate work.
To enforce good hygiene, i.e. to avoid adding apples and oranges, I think we should be restrictive on the quantities that can be added to/subtracted from Spectrum 1D's:
1) They should have matching flux units
2) They should have matching spectral axes and units.
For multiplication and division: 1) flux units may be present, or not. 2) They should have matching spectral axes and units.
While it may be an added burden to attach a spectral axis for the cases discussed by @rosteen above, it only takes one extra step.
@eteq @nmearl @keflavich Any thoughts you want to throw in the pile here?
Yeah, I agree with @pllim's suggestion that Spectrum1D
s should be required for non-scalars. This definitely makes it harder for users, so we need clear, explicit documentation about how to turn an array, or a model, into an appropriate Spectrum1D
.
Scalar addition & subtraction should be allowed with matching units. Scalar multiplication and division should be allowed with or without units.
Scalar multiplication and division should be allowed with or without units
What does it even mean physically when you multiply Jy by any random unit, like another Jy? Or Jy with Kelvin? Multiplication/division only makes sense to me when the other operand is unitless. Am I missing something?
What does it even mean physically when you multiply Jy by any random unit, like another Jy? Or Jy with Kelvin? Multiplication/division only makes sense to me when the other operand is unitless. Am I missing something?
You're right, that's a pretty rare case. There are some where I do that, though. For example, if you want to integrate over a spectrum, you might do so by multiplying by the spectral width of the pixel and then summing. If the pixels have varying width, you need to be able to multiply by a vector of pixel widths. So multiplying by, e.g., km/s or nm or Hz is sometimes reasonable.
There may be other cases - I don't think we should strictly exclude multiplication with unit'd quantities.
Ok, sounds like I'm outvoted on requiring Spectrum1D
for non-scalars 👍
I think I have implemented all the requests here. Please review https://github.com/astropy/specutils/pull/998 . The tests in #998 should give you a good idea of the expected new behaviors. Thanks!
Arithmetic operations succeed so long as the flux arrays are of compatible shape, but do not check for equivalent spectral axis values. For example,
succeeds but returns a
Spectrum1D
object with a spectral axis adopted directly fromsp2
. Instead this should raise an error (or perhaps eventually support optional interpolation for cases in which there is sufficient overlap).Related: #198