#35 Incorrect flux scaling

[sbourke]

Identified by mweastwood.

Fluxes appear to be scaled by sqrt(l^2 + m^2)
Test data with 1 Jy source at phase centre and 1 Jy source at 60 deg from phase centre:
http://rez.caltech.edu/~stephen/t1.ms.tgz
wsclean -size 512 512 -scale 0.25 t1.ms

Second source imaged at ~ 0.5 flux in wsclean. Maybe this is why PSFs fade out too?

[André Offringa]

Hi again sbourke ;-)

I've downloaded your set and imaged it, and can confirm that the WSClean images a 1 Jy and a 0.5 Jy source. However, when I overwrite your measurement set with the model below using an analytical evaluation of the visibility function, with the model containing two 100 Jy sources at the same positions as in your MS, I get two perfect 100 Jy sources out in WSClean. Hence, it seems the sources in your measurement set have not been predicted properly. In particular, I predict using Eq (1) from the WSClean paper, which comes from the radio book by Thompson et al. In that formula, there's a sqrt(1-l^2-m^2) term:

V = (A L) / sqrt(1-l^2-m^2) * exp(-2pi(uv + ...)..)

Maybe you have incorrectly predicted the visibilities without the sqrt(1-l^2-m^2) term...? The RMS of the PSF should go down away from the phase centre, so that's not necessarily wrong.

This is my sky model:

skymodel fileformat 1.1
source {
  name "Central"
  component {
    type point
    position 20h51m15s 36d29m40s
    measurement {
      frequency 74 MHz
      fluxdensity Jy 100.0 0 0 0.0
    }
  }
}
source {
  name "Northern"
  component {
    type point
    position 08h55m00s 80d39m25s
    measurement {
      frequency 74 MHz
      fluxdensity Jy 100.0 0 0 0.0
    }
  }
}

[sbourke]

Hi, hopefully I don't end up in your spam folder :-)

The sqrt(1-l^2-m^2) is a result of going from a spherical integral dΩ to dl dm and accounts for the sine projection. You're effectively applying it twice. In the analytical evaluation you're incorporating the compensating amplification in to the visibilities. If you rotate the phase centre all the fluxes will be wrong.

[André Offringa]

Hmmm I think you are right O_o.

I guess it's a bit complicated because the beam is so much different at high angles that flux levels in the dirty image are hard to interpret. The restoring beam should be smaller at high angles, but WSClean (as well as the other imagers I think) doesn't do that, and hence to keep the "integrated flux over constant solid angle" correctly it should apply the sqrt(1-l^2-m^2) term, and sources at high angle should have lower peak flux. However, on the other hand I notice that if you measure integrated flux with e.g. kvis it also doesn't take into account that the pixel size is larger at high angles, i.e. it calculates the integrated flux as if the pixel angle remains constant. The same probably holds for most source detectors (?). Unless of course the image is regridded...

I'll have to think & compare a bit what would be the best approach; I thought I'm using the most common way to handle this and I thought I had compared this with CASA. However for some reason CASA doesn't like the MS. I'll try another MS.

I you want to remove this normalization, you can do so with changes these two lines; in layeredimager.cpp line 535:

*dataPtr *= multiplicationFactor * sqrt(1.0 - l*l - m*m);

Should become:

*dataPtr *= multiplicationFactor;

And to also remove this inverse normalization during prediction (for Cotton-Schwab), in layeredimager.cpp line 595:

*dataPtr = *inPtr * multiplicationFactor / sqrt(1.0 - l*l - m*m);

Should become:

*dataPtr = *inPtr * multiplicationFactor;

[Michael Eastwood]

Hi Andre,

I guess it's a bit complicated because the beam is so much different at high angles that flux levels in the dirty image are hard to interpret. The restoring beam should be smaller at high angles, but WSClean (as well as the other imagers I think) doesn't do that, and hence to keep the "integrated flux over constant solid angle" correctly it should apply the sqrt(1-l^2-m^2) term, and sources at high angle should have lower peak flux. However, on the other hand I notice that if you measure integrated flux with e.g. kvis it also doesn't take into account that the pixel size is larger at high angles, i.e. it calculates the integrated flux as if the pixel angle remains constant. The same probably holds for most source detectors (?). Unless of course the image is regridded...

It's not really clear what you're trying to say here. Each pixel of the image should give an estimate of the flux per beam size in that direction. The peak flux of a 1 Jy point source should be 1 Jy / beam in the dirty image even as you move the source towards the horizon.

-- Michael

[André Offringa]

Hi Michael / sbourke,

I tested CASA's behaviour and it does not take out the sqrt(1-l^2-m^2) factor, so I'll remove it from WSClean too. As you say, that makes most sense, so that point sources have a peak flux which is correct.

My concern was keeping integrated flux values consistent too (...with other tools). [edit]

In any case, I agree that the sqrt(1-l^2-m^2) correction as it was does not make sense, so I have removed it. It will be fixed in WSClean 1.8. If you want the fix earlier: it's also already fixed on SF Git so you can do a git pull, or you can apply the changes I've detailed in my previous post.

Thanks for the report!