Open jballet opened 6 months ago
Thanks @jballet,
If I find what's causing this issue, how should I go about validating the correct values? Is the WCS representation (or are any other representations that I'm less aware of) also incorrect or is it just HEALPix?
Hi Alex, I have no reason to believe that the WCS version of gtsrcmaps (whatever projection) is incorrect. Everything was always perfectly consistent in that framework, only the HEALPix representation is at fault. So yes, the WCS projection that I included in the test harness can be considered the reference.
I've run into this issue too now, and I noticed, at least, that the scaling is different depending on whether it's an analytic model (RadialGaussian/RadialDisk) vs. a template model. If it's helpful, it seems like the Radial models are underscaled by a factor of about 896 (7128) or 900 (30^2) in case the numerolgy rings any bells. If I have time I'll try to figure out if the factor depends on the nside used in the binning.
Jean's test harness
I found a similar issue about "fermitools + healpix" a few years ago, and concluded that intensity by gtlike/gtmodel is given per deg2 instead of per sr, making a scale_factor larger by a factor of (180/pi)^2~3300. I have attached a note (in pdf) for reference. I guess the root cause is common, and hope this memo helps to interpret the scaling you found. Normalization issue in fermitools with Healpix map.pdf
I did some simulations to try to diagnose this a little more thoroughly. I went through the full toolchain to produce gtmodel maps for some extended sources using CAR projection (20x20 deg, 0.1 deg pixels) and healpix maps at order 5, 6, 7, 8, and 9. I found
This is a smoking gun that the problem is under resolving the convolution/sampling. What the other factor of 800ish is, I've got no ideas...
Here are the numerical results: 4FGL J0519.9-6845e: this is the LMC, a large map
4FGL J1552.4-5612e: this is MSH something or other, a small map
4FGL J1409.1-6121e: this is a RadialDisk source
4FGL J1834.5-0846e: this is a RadialGaussian source
Here is one more clue: I ran the same procedure but with refactor=8 in gtmodel. The default is 2, so I guess this means the internal pixels are 16x smaller in area. Lo and behold the offset here is smaller by a factor of 16. The improved consistency between the very coarse and the very fine HPX maps suggest to me this also helps the undersampling problem.
4FGL J0519.9-6845e CAR/HPX5=8.18 CAR/HPX6=32.72 CAR/HPX7=51.12 CAR/HPX8=51.12 CAR/HPX9=51.12
4FGL J1552.4-5612e CAR/HPX5=9.13 CAR/HPX6=37.42 CAR/HPX7=49.38 CAR/HPX8=50.36 CAR/HPX9=52.22
4FGL J1409.1-6121e CAR/HPX5=53.14 CAR/HPX6=55.39 CAR/HPX7=57.58 CAR/HPX8=57.60 CAR/HPX9=57.60
4FGL J1834.5-0846e CAR/HPX5=51.03 CAR/HPX6=52.22 CAR/HPX7=52.41 CAR/HPX8=52.61 CAR/HPX9=52.63
Jean-Marc Casandjian cannot enter extended sources into his all-sky modeling of diffuse emission, because gtsrcmaps returns a minuscule normalization (about 1000 times smaller than it should be) for extended sources in the HEALPix representation. I have put a test harness at this place. The file should be there until the end of 12/2023. The gtsrcmaps.sh script calls gtsrcmaps for one extended source and one nearby point source, in HEALPix then in a local (ARC) projection. The output files are there too. The maximum value in HEALPix and ARC is similar for the point source, but about 1000 times smaller in HEALPix than ARC for the extended source.