C++ Surface Brightness Profiles

SBProfile Base Class

class SBProfile

A base class representing all of the 2D surface brightness profiles that we know how to draw.

The SBProfile class is a representation of a surface brightness distribution across a 2-dimensional image plane, with real and/or Fourier-domain models of a wide variety of galaxy shapes, point-spread functions (PSFs), and their convolutions. There are several realizations of the SBProfile classes: There are the “atomic” classes that represent specific analytic profiles: (SBDeltaFunction, SBGaussian, SBSersic, SBAiry, SBExponential, SBBox, SBDeVaucouleurs and SBMoffat). SBInterpolatedImage represents a pattern defined by a grid of pixel values and a chosen interpolation scheme between pixel centers. SBTransform represents any affine transformation (shear, magnification, rotation, translation, and/or flux rescaling) of any other SBProfile. SBAdd represents the sum of any number of SBProfiles. SBConvolve represents the convolution of any number of SBProfiles, and SBDeconvolve is the deconvolution of one SBProfile with another.

Every SBProfile knows how to draw an Image<float> of itself in real and k space. Each also knows what is needed to prevent aliasing or truncation of itself when drawn. The classes model the surface brightness of the object, so the xValue routine returns the value of the surface brightness at a given point. Usually we define this as flux/arcsec^2.

However, when drawing an SBProfile, any distance measures used to define the SBProfile will be taken to be in units of pixels. Thus the image, which is nominally a surface brightness image is also a flux image. The flux in each pixel is the flux/pixel^2. The python layer takes care of the typical case of defining an SBProfile in terms of arcsec, then converting the image to pixels (scaling the size and flux appropriately), so that the draw command here has the correct flux.

This isn’t an abstract base class. An SBProfile is a concrete object which internally has a pointer to the implementation details (which is an abstract base class). Furthermore, all SBProfiles are immutable objects. Any changes are made through modifiers that return a new object. (e.g. setFlux, shear, shift, etc.) This means that we can safely make SBProfiles use shallow copies, since that will never be confusing, which in turn means that SBProfiles can be safely returned by value, used in containers (e.g. list<SBProfile>), etc.

The only constructor for SBProfile is the copy constructor. All SBProfiles need to be created as one of the derived types that have real constructors.

Well, technically, there is also a default constructor to make it easier to use containers of SBProfiles. However, it is an error to use an SBProfile that has been default constructed for any purpose.

The assignment operator does a shallow copy, replacing the current contents of the SBProfile with that of the rhs profile.

Subclassed by galsim::SBAdd, galsim::SBAiry, galsim::SBAutoConvolve, galsim::SBAutoCorrelate, galsim::SBBox, galsim::SBConvolve, galsim::SBDeconvolve, galsim::SBDeltaFunction, galsim::SBExponential, galsim::SBFourierSqrt, galsim::SBGaussian, galsim::SBInclinedExponential, galsim::SBInclinedSersic, galsim::SBInterpolatedImage, galsim::SBInterpolatedKImage, galsim::SBKolmogorov, galsim::SBMoffat, galsim::SBSecondKick, galsim::SBSersic, galsim::SBShapelet, galsim::SBSpergel, galsim::SBTopHat, galsim::SBTransform, galsim::SBVonKarman

Public Functions

SBProfile()

Default constructor for convenience only. Do not use!

This constructor is only provided so you can do things like:

std::list<SBProfile> prof_list;
prof_list.push_back(psf);
prof_list.push_back(gal);
prof_list.push_back(pix);
The default constructor for std::list strangely requires a default constructor for the argument type, even though it isn’t ever really used.

SBProfile(const SBProfile &rhs)

Only legitimate public constructor is a copy constructor.

SBProfile &operator=(const SBProfile &rhs)

operator= replaces the current contents with those of the rhs.

~SBProfile()

Destructor isn’t virtual, since derived classes don’t have anything to cleanup.

GSParams getGSParams() const

Get the GSParams object for this SBProfile Makes a copy, since this is used for python pickling, so the current SBProfile may go out of scope before we need to use the returned GSParams.

double xValue(const Position<double> &p) const

Return value of SBProfile at a chosen 2D position in real space.

Assume all are real-valued. xValue() may not be implemented for derived classes (SBConvolve) that require an FFT to determine real-space values. In this case, an SBError will be thrown.

Parameters:

p[in] 2D position in real space.

std::complex<double> kValue(const Position<double> &k) const

Return value of SBProfile at a chosen 2D position in k space.

Parameters:

k[in] 2D position in k space.

void getXRange(double &xmin, double &xmax, std::vector<double> &splits) const

Define the range over which the profile is not trivially zero.

These values are used when a real-space convolution is requested to define the appropriate range of integration. The implementation here is +- infinity for both x and y. Derived classes may override this if they a have different range.

void getYRange(double &ymin, double &ymax, std::vector<double> &splits) const
void getYRangeX(double x, double &ymin, double &ymax, std::vector<double> &splits) const
double maxK() const

Value of k beyond which aliasing can be neglected.

inline double nyquistDx() const

Real-space image pixel spacing that does not alias maxK.

double stepK() const

Sampling in k-space necessary to avoid folding too much of image in x space.

int getGoodImageSize(double dx) const

Determine a good size for a drawn image based on dx and stepK()

The basic formula is 2pi / (dx * stepK()). But then we round up to the next even integer value.

Parameters:

dx[in] The pixel scale of the image

Returns:

the recommended image size.

bool isAxisymmetric() const

Check whether the SBProfile is known to have rotational symmetry about x=y=0.

If the SBProfile has rotational symmetry, certain calculations can be simplified.

bool hasHardEdges() const

The presence of hard edges help determine whether real space convolution might be a better choice.

bool isAnalyticX() const

Check whether the SBProfile is analytic in the real domain.

An SBProfile is “analytic” in the real domain if values can be determined immediately at any position through formula or a stored table (no DFT); this makes certain calculations more efficient.

bool isAnalyticK() const

Check whether the SBProfile is analytic in the Fourier domain.

An SBProfile is “analytic” in the k domain if values can be determined immediately at any position through formula or a stored table (no DFT); this makes certain calculations more efficient.

Position<double> centroid() const

Returns (X, Y) centroid of SBProfile.

double getFlux() const

Get the total flux of the SBProfile.

double maxSB() const

Get an estimate of the maximum surface brightness.

SBTransform scaleFlux(double fluxRatio) const

Multiply the flux by fluxRatio.

SBTransform expand(double scale) const

Apply an overall scale change to the profile, preserving surface brightness.

SBTransform rotate(double theta) const

Apply a given rotation in radians.

SBTransform transform(double dudx, double dudy, double dvdx, double dvdy) const

Apply a transformation given by an arbitrary Jacobian matrix.

SBTransform shift(const Position<double> &delta) const

Apply a translation.

void shoot(PhotonArray &photons, BaseDeviate rng) const

Shoot photons through this SBProfile.

Returns an array of photon coordinates and fluxes that are drawn from the light distribution of this SBProfile. Absolute value of each photons’ flux should be approximately equal, but some photons can be negative as needed to represent negative regions. Note that the ray-shooting method is not intended to produce a randomized value of the total object flux, so do not assume that there will be sqrt(N) error on the flux. In fact most implementations will return a PhotonArray with exactly correct flux, with only the distribution of flux on the sky that will definitely have sampling noise.

The one definitive gaurantee is that, in the limit of large number of photons, the surface brightness distribution of the photons will converge on the SB pattern defined by the object.

Objects with regions of negative flux will result in creation of photons with negative flux. Absolute value of negative photons’ flux should be nearly equal to the standard flux of positive photons. Shot-noise fluctuations between the number of positive and negative photons will produce noise in the total net flux carried by the output PhotonArray.

The typical implementation will be to take the integral of the absolute value of flux, and divide it nearly equally into N photons. The photons are then drawn from the distribution of the absolute value of flux. If a photon is drawn from a region of negative flux, then that photon’s flux is negated. Because of cancellation, this means that each photon will carry more than getFlux()/N flux if there are negative-flux regions in the object. It also means that during convolution, addition, or interpolation, positive- and negative-flux photons can be contributing to the same region of the image. Their cancellation means that the shot noise may be substantially higher than you would expect if you had only positive-flux photons.

The photon flux may also vary slightly as a means of speeding up photon-shooting, as an alternative to rejection sampling. See OneDimensionalDeviate documentation.

Parameters:
  • photons[in] PhotonArray in which to write the photon information

  • rng[in] BaseDeviate that will be used to draw photons from distribution.

double getPositiveFlux() const

Return expectation value of flux in positive photons when shoot() is called.

Returns expectation value of flux returned in positive-valued photons when shoot() is called for this object. Default implementation is to return getFlux(), if it is positive, or 0 otherwise, which will be the case when the SBProfile is constructed entirely from elements of the same sign.

It should be generally true that getPositiveFlux() - getNegativeFlux() returns the same thing as getFlux(). Small difference may accrue from finite numerical accuracy in cases involving lookup tables, etc.

Returns:

Expected positive-photon flux.

double getNegativeFlux() const

Return expectation value of absolute value of flux in negative photons from shoot()

Returns expectation value of (absolute value of) flux returned in negative-valued photons when shoot() is called for this object. Default implementation is to return getFlux() if it is negative, 0 otherwise, which will be the case when the SBProfile is constructed entirely from elements that have the same sign.

It should be generally true that getPositiveFlux() - getNegativeFlux() returns the same thing as getFlux(). Small difference may accrue from finite numerical accuracy in cases involving lookup tables, etc.

Returns:

Expected absolute value of negative-photon flux.

template<typename T>
void draw(ImageView<T> image, double dx, double *jac, double xoff, double yoff, double flux_ratio) const

Draw an image of the SBProfile in real space forcing the use of real methods where we have a formula for x values. For SBProfiles without an analytic real-space representation, an exception will be thrown.

The image is not cleared out before drawing. So this profile will be added to anything already on the input image.

Parameters:
  • image[inout] (any of ImageViewF, ImageViewD, ImageViewS, ImageViewI)

  • dx, the[in] pixel scale

template<typename T>
void drawK(ImageView<std::complex<T>> image, double dk, double *jac) const

Draw an image of the SBProfile in k space.

For drawing in k space: routines are analagous to real space, except the image is complex. The image is normalized such that I(0,0) is the total flux.

Parameters:
  • image[inout] in k space (must be an ImageViewC)

  • dk, the[in] step size in k space

Simple C++ Profiles

class SBGaussian : public galsim::SBProfile

Gaussian Surface Brightness Profile.

The Gaussian Surface Brightness Profile is characterized by two properties, its flux and the characteristic size sigma where the radial profile of the circular Gaussian drops off as exp[-r^2 / (2. * sigma^2)]. The maxK() and stepK() are for the SBGaussian are chosen to extend to 4 sigma in both real and k domains, or more if needed to reach the folding_threshold spec.

Public Functions

SBGaussian(double sigma, double flux, const GSParams &gsparams)

Constructor.

Parameters:
  • sigma[in] characteristic size, surface brightness scales as exp[-r^2 / (2. * sigma^2)].

  • flux[in] flux of the Surface Brightness Profile.

  • gsparams[in] GSParams object storing constants that control the accuracy of image operations and rendering, if different from the default.

SBGaussian(const SBGaussian &rhs)

Copy constructor.

~SBGaussian()

Destructor.

double getSigma() const

Returns the characteristic size sigma of the Gaussian profile.

class SBDeltaFunction : public galsim::SBProfile

Delta Function Surface Brightness Profile.

The Delta Function Surface Brightness Profile is characterized by a single property, its ‘flux’.

Note that the DeltaFunction SBP cannot be drawn by itself. Instead, it should be applied as part of a convolution first.

Public Functions

SBDeltaFunction(double flux, const GSParams &gsparams)

Constructor.

Parameters:
  • flux[in] flux of the Surface Brightness Profile.

  • gsparams[in] GSParams object storing constants that control the accuracy of image operations and rendering, if different from the default.

SBDeltaFunction(const SBDeltaFunction &rhs)

Copy constructor.

~SBDeltaFunction()

Destructor.

class SBBox : public galsim::SBProfile

Surface Brightness Profile for the Boxcar function.

The boxcar function is a rectangular box. Convolution with a Boxcar function of dimensions width x height and sampling at pixel centres is equivalent to pixelation (i.e. Surface Brightness integration) across rectangular pixels of the same dimensions. This class is therefore useful for pixelating SBProfiles.

Public Functions

SBBox(double width, double height, double flux, const GSParams &gsparams)

Constructor.

Parameters:
  • width[in] width of Boxcar function along x.

  • height[in] height of Boxcar function along y.

  • flux[in] flux.

  • gsparams[in] GSParams object storing constants that control the accuracy of image operations and rendering, if different from the default.

SBBox(const SBBox &rhs)

Copy constructor.

~SBBox()

Destructor.

double getWidth() const

Returns the x dimension width of the Boxcar.

double getHeight() const

Returns the y dimension width of the Boxcar.

class SBTopHat : public galsim::SBProfile

Surface Brightness Profile for the TopHat function.

The tophat function is much like the boxcar, but a circular plateau, rather than a rectangle. It is defined by a radius and a flux.

Public Functions

SBTopHat(double radius, double flux, const GSParams &gsparams)

Constructor.

Parameters:
  • radius[in] radius of TopHat function

  • flux[in] flux.

  • gsparams[in] GSParams object storing constants that control the accuracy of image operations and rendering, if different from the default.

SBTopHat(const SBTopHat &rhs)

Copy constructor.

~SBTopHat()

Destructor.

double getRadius() const

Returns the radius of the TopHat.

PSF C++ Profiles

class SBMoffat : public galsim::SBProfile

Surface Brightness for the Moffat Profile (an approximate description of ground-based PSFs).

The Moffat surface brightness profile is I(R) propto [1 + (r/r_scale)^2]^(-beta). The SBProfile representation of a Moffat profile also includes an optional truncation beyond a given radius.

Public Functions

SBMoffat(double beta, double scale_radius, double trunc, double flux, const GSParams &gsparams)

Constructor.

Parameters:
  • beta[in] Moffat beta parameter for profile [1 + (r / rD)^2]^beta.

  • scale_radius[in] Scale radius, rD.

  • trunc[in] Outer truncation radius in same physical units as size, trunc = 0. for no truncation.

  • flux[in] Flux.

  • gsparams[in] GSParams object storing constants that control the accuracy of image operations and rendering, if different from the default.

SBMoffat(const SBMoffat &rhs)

Copy constructor.

~SBMoffat()

Destructor.

double getBeta() const

Returns beta of the Moffat profile [1 + (r / rD)^2]^beta.

double getFWHM() const

Returns the FWHM of the Moffat profile.

double getScaleRadius() const

Returns the scale radius rD of the Moffat profile [1 + (r / rD)^2]^beta.

double getHalfLightRadius() const

Returns the half-light radius of the Moffat profile.

double getTrunc() const

Returns the truncation radius.

class SBAiry : public galsim::SBProfile

Surface Brightness Profile for the Airy disk (perfect diffraction-limited PSF for a circular aperture), with central obscuration.

maxK() is set at the hard limit for Airy disks, stepK() makes transforms go to at least 5 lam/D or EE>(1-folding_threshold). Schroeder (10.1.18) gives limit of EE at large radius. This stepK could probably be relaxed, it makes overly accurate FFTs. Note x & y are in units of lambda/D here. Integral over area will give unity in this normalization.

Public Functions

SBAiry(double lam_over_D, double obscuration, double flux, const GSParams &gsparams)

Constructor.

Parameters:
  • lam_over_D[in] lam_over_D = (lambda * focal length) / (telescope diam) if arg is focal plane position, else lam_over_D = lambda / (telescope diam) if arg is in radians of field angle.

  • obscuration[in] linear dimension of central obscuration as fraction of pupil dimension.

  • flux[in] flux.

  • gsparams[in] GSParams object storing constants that control the accuracy of image operations and rendering, if different from the default.

SBAiry(const SBAiry &rhs)

Copy constructor.

~SBAiry()

Destructor.

double getLamOverD() const

Returns lam_over_D param of the SBAiry.

double getObscuration() const

Returns obscuration param of the SBAiry.

class SBKolmogorov : public galsim::SBProfile

Surface Brightness Profile for a Kolmogorov turbulent spectrum.

Public Functions

SBKolmogorov(double lam_over_r0, double flux, const GSParams &gsparams)

Constructor.

Parameters:
  • lam_over_r0[in] lambda / r0 in the physical units adopted (user responsible for consistency), where r0 is the Fried parameter. The FWHM of the Kolmogorov PSF is ~0.976 lambda/r0 (e.g., Racine 1996, PASP 699, 108). Typical values for the Fried parameter are on the order of 10 cm for most observatories and up to 20 cm for excellent sites. The values are usually quoted at lambda = 500 nm and r0 depends on wavelength as [r0 ~ lambda^(6/5)].

  • flux[in] Flux.

  • gsparams[in] GSParams object storing constants that control the accuracy of image operations and rendering, if different from the default.

SBKolmogorov(const SBKolmogorov &rhs)

Copy constructor.

~SBKolmogorov()

Destructor.

double getLamOverR0() const

Returns lam_over_r0 param of the SBKolmogorov.

class SBVonKarman : public galsim::SBProfile

Public Functions

SBVonKarman(double lam, double r0, double L0, double flux, double scale, bool doDelta, const GSParams &gsparams, double force_stepk)

Constructor.

Parameters:
  • lam[in] Wavelength in nm.

  • r0[in] Fried parameter in m (at given wavelength lam).

  • L0[in] Outer scale in m.

  • flux[in] Flux.

  • scale[in] Scale of ‘x’ in xValue in arcsec.

  • doDelta[in] Whether or not to include delta-function contribution to encircled energy when computing stepk/maxk/HLR.

  • gsparams[in] GSParams.

SBVonKarman(const SBVonKarman &rhs)

Copy constructor.

~SBVonKarman()

Destructor.

double getLam() const

Getters.

double getR0() const
double getL0() const
double getScale() const
bool getDoDelta() const
double getDelta() const
double getHalfLightRadius() const
double structureFunction(double) const

Friends

friend class VKXIntegrand
class SBSecondKick : public galsim::SBProfile

Public Functions

SBSecondKick(double lam_over_r0, double kcrit, double flux, const GSParamsPtr &gsparams)

Constructor.

Parameters:
  • lam_over_r0[in] lambda/r0, equivalent to the same in SBKolmogorov

  • kcrit[in] Critical turbulence Fourier mode in units of r0.

  • flux[in] Flux.

  • gsparams[in] GSParams.

SBSecondKick(const SBSecondKick &rhs)

Copy constructor.

~SBSecondKick()

Destructor.

double getLamOverR0() const

Getters.

double getKCrit() const
double getDelta() const
double kValue(double) const

Alternate versions of x/k Value for testing purposes.

double kValueRaw(double) const
double xValue(double) const
double xValueRaw(double) const
double xValueExact(double) const
double structureFunction(double) const

Galaxy C++ Profiles

class SBExponential : public galsim::SBProfile

Exponential Surface Brightness Profile.

Surface brightness profile with I(r) propto exp[-r/r_0] for some scale-length r_0. This is a special case of the Sersic profile, but is given a separate class since the Fourier transform has closed form and can be generated without lookup tables.

Public Functions

SBExponential(double r0, double flux, const GSParams &gsparams)

Constructor - note that r0 is scale length, NOT half-light radius re as in SBSersic.

Parameters:
  • r0[in] scale length for the profile that scales as exp[-(r / r0)], NOT the half-light radius re.

  • flux[in] flux.

  • gsparams[in] GSParams object storing constants that control the accuracy of image operations and rendering, if different from the default.

SBExponential(const SBExponential &rhs)

Copy constructor.

~SBExponential()

Destructor.

double getScaleRadius() const

Returns the scale radius of the Exponential profile.

class SBSersic : public galsim::SBProfile

Sersic Surface Brightness Profile.

The Sersic Surface Brightness Profile is characterized by three properties: its Sersic index n, its flux, and the half-light radius re (or scale radius r0). Given these properties, the surface brightness profile scales as I(r) propto exp[-(r/r0)^{1/n}], or I(r) propto exp[-b*(r/re)^{1/n}]. The code is limited to 0.3 <= n <= 6.2, with an exception thrown for values outside that range.

The SBProfile representation of a Sersic profile also includes an optional truncation beyond a given radius, by the parameter trunc. Internal Sersic information are cached according to the (n, trunc/r0) pair. All internal calculations are based on the scale radius r0. If the profile is specified by the half-light radius re, the corresponding scale radius r0 is calculated, and vice versa.

When the Sersic profile is specfied by the scale radius with truncation, the normalization is adjusted such that the truncated profile has the specified flux (its half-light radius will differ from an equivalent Sersic without truncation). Similarly, when the Sersic profile is specified by the half-light radius with truncation, SBSersic generates a profile whose flux and half-light radius is as specified, by adjusting its normalization and scale radius.

Another optional parameter, flux_untruncated = true, allows the setting of the flux to the untruncated Sersic, while generating a truncated Sersic (i.e., the normalizaton is the same with respect to the untruncated case). This facilitates the comparison of truncated and untruncated Sersic, as the amplitude (as well as the scale radius r0, if half-light radius is specified) changes when a truncated Sersic is specified with flux_untruncated = false. The flux_untruncated variable is ignored if trunc = 0.

Note that when trunc > 0. and flux_untruncated == true, the actual flux will not be the same as the specified value; its true flux is returned by the getFlux() method. Similarly for the half-light radius, when the Sersic profile is specified by the half-light radius; the getHalfLightRadius() method will return the true half-light radius. The scale radius will remain at the same value, if this quantity was used to specify the profile.

There are two special cases of the Sersic profile that have their own SBProfiles: n=1 (SBExponential), n=0.5 (SBGaussian). These special cases use several simplifications in all calculations, whereas for general n, the Fourier transform must be treated numerically.

Public Functions

SBSersic(double n, double scale_radius, double flux, double trunc, const GSParams &gsparams)

Constructor.

Parameters:
  • n[in] Sersic index.

  • scale_radius[in] Scale radius

  • flux[in] Flux.

  • trunc[in] Outer truncation radius in same physical units as size; trunc = 0. for no truncation.

  • gsparams[in] GSParams object storing constants that control the accuracy of image operations and rendering, if different from the default.

SBSersic(const SBSersic &rhs)

Copy constructor.

~SBSersic()

Destructor.

double getN() const

Returns the Sersic index n of the profile.

double getScaleRadius() const

Returns the scale radius r0 of the Sersic profile exp[-(r/r_0)^(1/n)].

double getHalfLightRadius() const

Returns the half light radius of the Sersic profile.

double getTrunc() const

Returns the truncation radius.

class SBInclinedExponential : public galsim::SBProfile

Inclined exponential surface brightness profile.

Surface brightness profile based on the distribution I(R,z) propto sech^2(z/Hs)*exp(-R/Rs), where Hs is the scale height of the disk, Rs is the scale radius, z is the distance along the minor axis, and R is the distance perpendicular to it. This profile is determined by four parameters: The inclination angle, the scale radius, the scale height, and the flux.

Note that the position angle is always zero. A profile with a different position angle can be obtained through the rotate() method of the corresponding Python class.

Public Functions

SBInclinedExponential(double inclination, double scale_radius, double scale_height, double flux, const GSParams &gsparams)

Constructor.

Parameters:
  • inclination[in] Inclination angle, where 0 = face-on, pi/2 = edge-on

  • scale_radius[in] Scale radius of the exponential disk.

  • scale_height[in] Scale height of the exponential disk.

  • flux[in] Flux.

  • gsparams[in] GSParams object storing constants that control the accuracy of image operations and rendering, if different from the default.

SBInclinedExponential(const SBInclinedExponential &rhs)

Copy constructor.

~SBInclinedExponential()

Destructor.

double getInclination() const

Returns the inclination angle of the profile in radians.

double getScaleRadius() const

Returns the scale radius r0 of the disk profile.

double getScaleHeight() const

Returns the scale height h0 of the disk profile.

class SBInclinedSersic : public galsim::SBProfile

Inclined Sersic Surface Brightness Profile.

This class is similar to a Sersic profile, additionally allowing inclination relative to the line of sight. The true profile is assumed to follow a Sersic distribution in R, multiplied by sech^2(z/Hs), where Hs is the scale height of the disk and z is the distance along the minor axis. The inclination angle determines how elliptical the profile appears.

See the documentation for the SBSersic class for further details on Sersic profiles.

Note that the position angle is always zero. A profile with a different position angle can be obtained through the rotate() method of the corresponding Python class.

If the inclination will always be zero (face-on), the SBSersic class can instead be used as a slightly faster alternative. If no truncation radius will be applied and n=1, the SBInclinedExponential class can be used as a much faster alternative.

Public Functions

SBInclinedSersic(double n, double inclination, double scale_radius, double height, double flux, double trunc, const GSParams &gsparams)

Constructor.

Parameters:
  • n[in] Sersic index.

  • inclination[in] Inclination of the disk relative to line of sight, in radians, where 0 = face-on and pi/2 = edge-on.

  • scale_radius[in] Scale radius

  • height[in] Scale height

  • flux[in] Flux.

  • trunc[in] Outer truncation radius in same physical units as size; trunc = 0. for no truncation.

  • gsparams[in] GSParams object storing constants that control the accuracy of image operations and rendering, if different from the default.

SBInclinedSersic(const SBInclinedSersic &rhs)

Copy constructor.

~SBInclinedSersic()

Destructor.

double getN() const

Returns the Sersic index n of the profile.

double getInclination() const

Returns the inclination angle of the profile in radians.

double getScaleRadius() const

Returns the scale radius r0 of the Sersic profile.

double getHalfLightRadius() const

Returns the half light radius of the Sersic profile (if it were face-on).

double getScaleHeight() const

Returns the scale height h0 of the disk profile.

double getTrunc() const

Returns the truncation radius.

class SBSpergel : public galsim::SBProfile

Spergel Surface Brightness Profile.

Surface brightness profile with I(r) propto (r/r_c)^{nu} K_{nu}(r/r_c) for some scale-length r_c = r_0 / c_{nu}, where K_{nu}(u) is the modified Bessel function of the second kind (also confusingly referred to as the ‘spherical modified Bessel function of the third kind’) and nu > -1. For different parameters nu this profile can approximate Sersic profiles with different indices.

Reference: D. N. Spergel, “ANALYTICAL GALAXY PROFILES FOR PHOTOMETRIC AND LENSING ANALYSIS,”” ASTROPHYS J SUPPL S 191(1), 58-65 (2010) [doi:10.1088/0067-0049/191/1/58].

Public Functions

SBSpergel(double nu, double scale_radius, double flux, const GSParams &gsparams)

Constructor - note that r0 is scale length, NOT half-light radius re as in SBSersic.

Parameters:
  • nu[in] index parameter setting the logarithmic slope of the profile.

  • scale_radius[in] scale radius

  • flux[in] flux.

  • gsparams[in] GSParams object storing constants that control the accuracy of image operations and rendering, if different from the default.

SBSpergel(const SBSpergel &rhs)

Copy constructor.

~SBSpergel()

Destructor.

double getNu() const

Returns the Spergel index nu of the profile.

double getScaleRadius() const

Returns the scale radius of the Spergel profile.

double calculateIntegratedFlux(double r) const

Return integrated flux of circular profile.

double calculateFluxRadius(double f) const

Return radius enclosing flux f.

Arbitrary C++ Profiles

class SBInterpolatedImage : public galsim::SBProfile

Surface Brightness Profile represented by interpolation over one or more data tables/images.

The SBInterpolatedImage class represents an arbitrary surface brightness profile (supplied as an image), including rules for how to interpolate the profile between the supplied pixel values.

It is assumed that input images oversample the profiles they represent. maxK() is set at the Nyquist frequency of the input image, although it should be noted that interpolants other than the ideal sinc function may make the max frequency higher than this. The output is required to be periodic on a scale > original image extent + kernel footprint, and stepK() is set accordingly.

The normal way to make an SBInterpolatedImage is to provide the image to interpolate and the interpolation scheme. See Interpolant.h for more about the different kind of interpolation.

You can provide different interpolation schemes for real and Fourier space (passed as xInterp and kInterp respectively). These are required, but there are sensible defaults in the python layer wrapper class, InterpolatedImage.

The ideal k-space interpolant is a sinc function; however, the quintic interpolant is the default, based on detailed investigations on the tradeoffs between accuracy and speed. Note that, as in Bernstein & Gruen (2012), the accuracy achieved by this interpolant is dependent on our choice of 4x pad factor. Users who do not wish to pad the arrays to this degree may need to use a higher-order Lanczos interpolant instead, but this is not the recommended usage. (Note: this padding is done by the python layer now, not here.)

The surface brightness profile will be in terms of the image pixels. The python layer InterpolatedImage class takes care of converting between these units and the arcsec units that are usually desired.

Public Functions

SBInterpolatedImage(const BaseImage<double> &image, const Bounds<int> &init_bounds, const Bounds<int> &nonzero_bounds, const Interpolant &xInterp, const Interpolant &kInterp, double stepk, double maxk, const GSParams &gsparams)

Initialize internal quantities and allocate data tables based on a supplied 2D image.

Parameters:
  • image[in] Input Image (ImageF or ImageD).

  • init_bounds[in] The bounds of the original unpadded image.

  • nonzero_bounds[in] The bounds in which the padded image is non-zero.

  • xInterp[in] Interpolation scheme to adopt between pixels

  • kInterp[in] Interpolation scheme to adopt in k-space

  • stepk[in] If > 0, force stepk to this value.

  • maxk[in] If > 0, force maxk to this value.

  • gsparams[in] GSParams object storing constants that control the accuracy of image operations and rendering.

SBInterpolatedImage(const SBInterpolatedImage &rhs)

Copy Constructor.

~SBInterpolatedImage()

Destructor.

const Interpolant &getXInterp() const
const Interpolant &getKInterp() const
double getPadFactor() const
void calculateStepK(double max_stepk = 0.) const

Refine the value of stepK if the input image was larger than necessary.

Parameters:

max_stepk[in] Optional maximum value of stepk if you have some a priori knowledge about an appropriate maximum.

void calculateMaxK(double max_maxk = 0.) const

Refine the value of maxK if the input image had a smaller scale than necessary.

Parameters:

max_maxk[in] Optional maximum value of maxk if you have some a priori knowledge about an appropriate maximum.

ConstImageView<double> getPaddedImage() const
ConstImageView<double> getNonZeroImage() const
ConstImageView<double> getImage() const
class SBInterpolatedKImage : public galsim::SBProfile

Public Functions

SBInterpolatedKImage(const BaseImage<std::complex<double>> &kimage, double stepk, const Interpolant &kInterp, const GSParams &gsparams)

Initialize internal quantities and allocate data tables based on a supplied 2D image.

Parameters:
  • kimage[in] Input Fourier-space Image (ImageC).

  • stepk[in] If > 0, force stepk to this value.

  • kInterp[in] Interpolation scheme to adopt in k-space

  • gsparams[in] GSParams object storing constants that control the accuracy of image operations and rendering.

SBInterpolatedKImage(const BaseImage<double> &data, double stepk, double maxk, const Interpolant &kInterp, const GSParams &gsparams)
SBInterpolatedKImage(const SBInterpolatedKImage &rhs)

Copy Constructor.

~SBInterpolatedKImage()

Destructor.

const Interpolant &getKInterp() const
ConstImageView<double> getKData() const
class SBShapelet : public galsim::SBProfile

Class for describing polar shapelet surface brightness profiles.

Public Functions

SBShapelet(double sigma, LVector bvec, const GSParams &gsparams)

Constructor.

Parameters:
  • sigma[in] scale size of Gauss-Laguerre basis set.

  • bvec[in] bvec[n,m] contains flux information for the (n, m) basis function.

  • gsparams[in] GSParams object storing constants that control the accuracy of image operations and rendering, if different from the default.

SBShapelet(const SBShapelet &rhs)

Copy Constructor.

~SBShapelet()

Destructor.

double getSigma() const
const LVector &getBVec() const
void rotate(double theta)

Composite C++ Profiles

class SBAdd : public galsim::SBProfile

Sums SBProfiles.

The SBAdd class can be used to add arbitrary numbers of SBProfiles together.

Public Functions

SBAdd(const std::list<SBProfile> &slist, const GSParams &gsparams)

Constructor, list of inputs.

Parameters:
  • slist[in] List of SBProfiles.

  • gsparams[in] GSParams object storing constants that control the accuracy of image operations and rendering, if different from the default.

SBAdd(const SBAdd &rhs)

Copy constructor.

~SBAdd()

Destructor.

std::list<SBProfile> getObjs() const

Get the list of SBProfiles that are being added together.

class SBConvolve : public galsim::SBProfile

Convolve SBProfiles.

Convolve two, three or more SBProfiles together.

The profiles to be convolved may be provided either as the first 2 or 3 parameters in the constructor, or as a std::list<SBProfile>.

The convolution will normally be done using discrete Fourier transforms of each of the component profiles, multiplying them together, and then transforming back to real space. The nominal flux of the resulting SBConvolve is the product of the fluxes of each of the component profiles. Thus, when using the SBConvolve to convolve a galaxy of some desired flux with a PSF, it is important to normalize the flux in the PSF to 1 beforehand.

The stepK used for the k-space image will be (Sum 1/stepK()^2)^(-1/2) where the sum is over all the components being convolved. Since the size of the convolved image scales roughly as the quadrature sum of the components, this should be close to Pi/Rmax where Rmax is the radius that encloses all but (1-folding_threshold) of the flux in the final convolved image.

The maxK used for the k-space image will be the minimum of the maxK() calculated for each component. Since the k-space images are multiplied, if one of them is essentially zero beyond some k value, then that will be true of the final image as well.

There is also an option to do the convolution as integrals in real space. Each constructor has an optional boolean parameter, real_space, that comes immediately after the list of profiles to convolve. Currently, the real-space integration is only enabled for 2 profiles. If you try to use it for more than 2 profiles, an exception will be thrown.

The real-space convolution is normally slower than the DFT convolution. The exception is if both component profiles have hard edges (e.g. a truncated Moffat with a Box). In that case, the maxK for each component is quite large since the ringing dies off fairly slowly. So it can be quicker to use real-space convolution instead.

Public Functions

SBConvolve(const std::list<SBProfile> &slist, bool real_space, const GSParams &gsparams)

Constructor, list of inputs.

Parameters:
  • slist[in] Input: list of SBProfiles.

  • real_space[in] Do convolution in real space?

  • gsparams[in] GSParams object storing constants that control the accuracy of image operations and rendering, if different from the default.

SBConvolve(const SBConvolve &rhs)

Copy constructor.

~SBConvolve()

Destructor.

std::list<SBProfile> getObjs() const

Get the list of SBProfiles that are being convolved.

bool isRealSpace() const

Return whether the convolution should be done in real space.

class SBAutoConvolve : public galsim::SBProfile

Public Functions

SBAutoConvolve(const SBProfile &s, bool real_space, const GSParams &gsparams)

Constructor.

Parameters:
  • s[in] SBProfile to be convolved with itself.

  • real_space[in] Do convolution in real space?

  • gsparams[in] GSParams to use, if different from the default.

SBAutoConvolve(const SBAutoConvolve &rhs)

Copy constructor.

~SBAutoConvolve()

Destructor.

SBProfile getObj() const

Get the SBProfile being convolved.

bool isRealSpace() const

Return whether the convolution should be done in real space.

class SBAutoCorrelate : public galsim::SBProfile

Public Functions

SBAutoCorrelate(const SBProfile &s, bool real_space, const GSParams &gsparams)

Constructor.

Parameters:
  • s[in] SBProfile to be correlated with itself.

  • real_space[in] Do convolution in real space?

  • gsparams[in] GSParams to use, if different from the default.

SBAutoCorrelate(const SBAutoCorrelate &rhs)

Copy constructor.

~SBAutoCorrelate()

Destructor.

SBProfile getObj() const

Get the SBProfile being conrrelated.

bool isRealSpace() const

Return whether the convolution should be done in real space.

Transformed C++ Profiles

class SBTransform : public galsim::SBProfile

An affine transformation of another SBProfile.

Origin of original shape will now appear at _cen. Flux is NOT conserved in transformation - surface brightness is preserved. We keep track of all distortions in a 2x2 matrix M = [(A B), (C D)] = [row1, row2] plus a 2-element Positon object cen for the shift, and a flux scaling, in addition to the scaling implicit in the matrix M = abs(det(M)).

Public Functions

SBTransform(const SBProfile &sbin, const double *jac, const Position<double> &cen, double ampScaling, const GSParams &gsparams)

General constructor.

Parameters:
  • obj[in] SBProfile being transformed

  • jac[in] 4-element array (A,B,C,D) of 2x2 distortion matrix M = [(A B), (C D)]

  • cen[in] 2-element (x, y) Position for the translational shift.

  • ampScaling[in] Amount by which the SB amplitude should be multiplied.

  • gsparams[in] GSParams object storing constants that control the accuracy of image operations and rendering, if different from the default.

SBTransform(const SBTransform &rhs)

Copy constructor.

~SBTransform()

Destructor.

SBProfile getObj() const
void getJac(double &mA, double &mB, double &mC, double &mD) const
Position<double> getOffset() const
double getFluxScaling() const
class SBDeconvolve : public galsim::SBProfile

SBProfile adapter which inverts its subject in k space to effect a deconvolvution.

Param adaptee:

[in] SBProfile for which to effect a deconvolution by k space inversion.

Param gsparams:

[in] GSParams object storing constants that control the accuracy of image operations and rendering, if different from the default.

Public Functions

SBDeconvolve(const SBProfile &adaptee, const GSParams &gsparams)

Constructor.

SBDeconvolve(const SBDeconvolve &rhs)

Copy constructor.

~SBDeconvolve()

Destructor.

SBProfile getObj() const

Get the SBProfile being deconvolved.

class SBFourierSqrt : public galsim::SBProfile

SBProfile adapter which computes the square root of its subject in k space.

Param adaptee:

[in] SBProfile to compute the Fourier-space square root of.

Param gsparams:

[in] GSParams object storing constants that control the accuracy of image operations and rendering, if different from the default.

Public Functions

SBFourierSqrt(const SBProfile &adaptee, const GSParams &gsparams)

Constructor.

SBFourierSqrt(const SBFourierSqrt &rhs)

Copy constructor.

SBProfile getObj() const

Get the SBProfile being operated on.

~SBFourierSqrt()

Destructor.

C++ GSParams

struct GSParams

Public Functions

GSParams(int _minimum_fft_size, int _maximum_fft_size, double _folding_threshold, double _stepk_minimum_hlr, double _maxk_threshold, double _kvalue_accuracy, double _xvalue_accuracy, double _table_spacing, double _realspace_relerr, double _realspace_abserr, double _integration_relerr, double _integration_abserr, double _shoot_accuracy)

A set of numbers that govern how SBProfiles make various speed/accuracy tradeoff decisions.

These parameters can be broadly split into two groups: i) parameters that affect the rendering of objects by Discrete Fourier Transform (DFT) and by real space convolution; and ii) parameters that affect rendering by photon shooting.

The DFT and real space convolution relevant params are:

The Photon Shooting relevant params are:

Parameters:
  • minimum_fft_size – Constant giving minimum FFT size we’re willing to do.

  • maximum_fft_size – Constant giving maximum FFT size we’re willing to do.

  • folding_threshold – A threshold parameter used for setting the stepK value for FFTs. The FFT’s stepK is set so that at most a fraction folding_threshold of the flux of any profile is folded.

  • stepk_minimum_hlr – In addition to the above constraint for aliasing, also set stepk such that pi/stepk is at least stepk_minimum_hlr times the profile’s half-light radius (for profiles that have a well-defined half-light radius).

  • maxk_threshold – A threshold parameter used for setting the maxK value for FFTs. The FFT’s maxK is set so that the k-values that are excluded off the edge of the image are less than maxk_threshold.

  • kvalue_accuracy – Accuracy of values in k-space. If a k-value is less than kvalue_accuracy, then it may be set to zero. Similarly, if an alternate calculation has errors less than kvalue_accuracy, then it may be used instead of an exact calculation. Note: This does not necessarily imply that all kvalues are this accurate. There may be cases where other choices we have made lead to errors greater than this. But whenever we do an explicit calculation about this, this is the value we use. This should typically be set to a lower, more stringent value than maxk_threshold.

  • xvalue_accuracy – Accuracy of values in real space. If a value in real space is less than xvalue_accuracy, then it may be set to zero. Similarly, if an alternate calculation has errors less than xvalue_accuracy, then it may be used instead of an exact calculation.

  • table_spacing – Several profiles use lookup tables for either the Hankel transform (Sersic, truncated Moffat) or the real space radial function (Kolmogorov). We try to estimate a good spacing between values in the lookup tables based on either xvalue_accuracy or kvalue_accuracy as appropriate. However, you may change the spacing with table_spacing. Using table_spacing < 1 will use a spacing value that much smaller than the default, which should produce more accurate interpolations.

  • realspace_relerr – The target relative accuracy for real-space convolution.

  • realspace_abserr – The target absolute accuracy for real-space convolution.

  • integration_relerr – Target relative accuracy for integrals (other than real-space convolution).

  • integration_abserr – Target absolute accuracy for integrals (other than real-space convolution).

  • shoot_accuracy – Accuracy of total flux for photon shooting. The photon shooting algorithm sometimes needs to sample the radial profile out to some value. We choose the outer radius such that the integral encloses at least (1-shoot_accuracy) of the flux.

inline GSParams()

A reasonable set of default values

bool operator==(const GSParams &rhs) const
bool operator<(const GSParams &rhs) const

Public Members

int minimum_fft_size
int maximum_fft_size
double folding_threshold
double stepk_minimum_hlr
double maxk_threshold
double kvalue_accuracy
double xvalue_accuracy
double table_spacing
double realspace_relerr
double realspace_abserr
double integration_relerr
double integration_abserr
double shoot_accuracy