---

8.1.1.5 Sampling from a function ¶

A pixel is the ultimate level of accuracy to gather data, we cannot get any more accurate in one image, this is known as sampling in signal processing. However, the mathematical profiles which describe our models have infinite accuracy. Over a large fraction of the area of astrophysically interesting profiles (for example, galaxies or PSFs), the variation of the profile over the area of one pixel is not too significant. In such cases, the elliptical radius (\(r_{el}\)) of the center of the pixel can be assigned as the final value of the pixel, (see Defining an ellipse and ellipsoid).

As you approach their center, some galaxies become very sharp (their value significantly changes over one pixel’s area). This sharpness increases with smaller effective radius and larger Sérsic values. Thus rendering the central value extremely inaccurate. The first method that comes to mind for solving this problem is integration. The functional form of the profile can be integrated over the pixel area in a 2D integration process. However, unfortunately numerical integration techniques also have their limitations and when such sharp profiles are needed they can become extremely inaccurate.

The most accurate method of sampling a continuous profile on a discrete space is by choosing a large number of random points within the boundaries of the pixel and taking their average value (or Monte Carlo integration). This is also, generally speaking, what happens in practice with the photons on the pixel. The number of random points can be set with --numrandom.

Unfortunately, repeating this Monte Carlo process would be extremely time and CPU consuming if it is to be applied to every pixel. In order to not loose too much accuracy, in MakeProfiles, the profile is built using both methods explained below. The building of the profile begins from its central pixel and continues (radially) outwards. Monte Carlo integration is first applied (which yields \(F_r\)), then the central pixel value (\(F_c\)) is calculated on the same pixel. If the fractional difference (\(|F_r-F_c|/F_r\)) is lower than a given tolerance level (specified with --tolerance) MakeProfiles will stop using Monte Carlo integration and only use the central pixel value.

The ordering of the pixels in this inside-out construction is based on \(r=\sqrt{(i_c-i)^2+(j_c-j)^2}\), not \(r_{el}\), see Defining an ellipse and ellipsoid. When the axis ratios are large (near one) this is fine. But when they are small and the object is highly elliptical, it might seem more reasonable to follow \(r_{el}\) not \(r\). The problem is that the gradient is stronger in pixels with smaller \(r\) (and larger \(r_{el}\)) than those with smaller \(r_{el}\). In other words, the gradient is strongest along the minor axis. So if the next pixel is chosen based on \(r_{el}\), the tolerance level will be reached sooner and lots of pixels with large fractional differences will be missed.

Monte Carlo integration uses a random number of points. Thus, every time you run it, by default, you will get a different distribution of points to sample within the pixel. In the case of large profiles, this will result in a slight difference of the pixels which use Monte Carlo integration each time MakeProfiles is run. To have a deterministic result, you have to fix the random number generator properties which is used to build the random distribution. This can be done by setting the GSL_RNG_TYPE and GSL_RNG_SEED environment variables and calling MakeProfiles with the --envseed option. To learn more about the process of generating random numbers, see Generating random numbers.

The seed values are fixed for every profile: with --envseed, all the profiles have the same seed and without it, each will get a different seed using the system clock (which is accurate to within one microsecond). The same seed will be used to generate a random number for all the sub-pixel positions of all the profiles. So in the former, the sub-pixel points checked for all the pixels undergoing Monte carlo integration in all profiles will be identical. In other words, the sub-pixel points in the first (closest to the center) pixel of all the profiles will be identical with each other. All the second pixels studied for all the profiles will also receive an identical (different from the first pixel) set of sub-pixel points and so on. As long as the number of random points used is large enough or the profiles are not identical, this should not cause any systematic bias.