In this vignette we will demonstrate fitting in the most automated mode available in ProFit. This is very high level, but hopefully usefully so in many cases. Users who need finer level control should look at some of the other vignettes that break various phases we do automatically here into more manual and explicit steps.
Get the latest version of ProFound and ProFit:
library(devtools)
install_github('ICRAR/ProFit')
install_github('asgr/ProFound')
install_github('asgr/Highlander')
Set global evaluate (basically TRUE for GitHub version and FALSE for CRAN):
= TRUE evalglobal
First load the libraries we need:
library(ProFit)
#> Loading required package: FITSio
library(ProFound)
#> Loading required package: Rcpp
library(magicaxis)
library(Highlander)
Load the data:
= readFITS(system.file("extdata", 'VIKING/mystery_VIKING_Z.fits', package="ProFound"))$imDat image
And take a look at what we have got:
magimage(image, hicut=1)
It is Z-band data from the VIKING survey. You can see that there are two very bright stars next to our galaxy of interest. At least the nearer of the two is producing enough extended flux that you profiling might be compromised unless we account for it.
If we do not have a good model for the PSF, we can use the high level fairly automated function profitAllStarDoFit to find multiple reasonable stars within the frame (6 by default) and fit them all with a common Moffat model for our PSF. This takes about 7 minutes to run on a modern single CPU.
= profitAllStarDoFit(image, magzero=30, SBdilate=2, rough=TRUE)
stars_auto #> Running Found2Fit
#> Running ProFound
#> Running Highander
#> Iteration 1
#> CMA 1: -1390.597 25.472 77.021 26.352 76.413 25.722 25.327 76.394 77.211 20.875 21.897 22.154 22.366 0.619 0.391
#>
#> Laplace's Demon was called on Tue Feb 23 17:52:08 2021
#>
#> Performing initial checks...
#> WARNING: Matrix mask may be rank-deficient.
#> WARNING: Matrix segim may be rank-deficient.
#> WARNING: Matrix region may be rank-deficient.
#> WARNING: Matrix calcregion may be rank-deficient.
#> Algorithm: Componentwise Hit-And-Run Metropolis
#>
#> Laplace's Demon is beginning to update...
#> Iteration: 100, Proposal: Componentwise, LP: -1392.7
#>
#> Assessing Stationarity
#> Assessing Thinning and ESS
#> Creating Summaries
#> Creating Output
#>
#> Laplace's Demon has finished.
#> Iteration 2
#>
#> Laplace's Demon was called on Tue Feb 23 17:52:13 2021
#>
#> Performing initial checks...
#> WARNING: Matrix mask may be rank-deficient.
#> WARNING: Matrix segim may be rank-deficient.
#> WARNING: Matrix region may be rank-deficient.
#> WARNING: Matrix calcregion may be rank-deficient.
#> Algorithm: Componentwise Hit-And-Run Metropolis
#>
#> Laplace's Demon is beginning to update...
#> Iteration: 100, Proposal: Componentwise, LP: -1394.3
#>
#> Assessing Stationarity
#> Assessing Thinning and ESS
#> Creating Summaries
#> Creating Output
#>
#> Laplace's Demon has finished.
And check the resulting fit:
profitLikeModel(stars_auto$parm, stars_auto$Data, makeplots = TRUE, plotchisq = TRUE)
#> Warning in dt(x, tdof): NaNs produced
#> $LP
#> [1] -1390.597
#>
#> $Dev
#> [1] 2781.194
#>
#> $Monitor
#> LL LP tend.elapsed
#> -1390.597 -1390.597 12.577
#>
#> $yhat
#> [1] 1
#>
#> $parm
#> moffat.xcen1 moffat.xcen2 moffat.xcen3 moffat.xcen4 moffat.ycen1 moffat.ycen2
#> 25.4717075 77.0208849 26.3522732 76.4128888 25.7223373 25.3272761
#> moffat.ycen3 moffat.ycen4 moffat.mag1 moffat.mag2 moffat.mag3 moffat.mag4
#> 76.3944542 77.2109552 20.8751854 21.8968733 22.1541462 22.3662960
#> moffat.fwhm1 moffat.con1
#> 0.6191344 0.3910980
And we can look at the model PSF too:
magimage(stars_auto$psf)
And look at the model parameters too (if you are interested):
$psf_modellist
stars_auto#> $moffat
#> $moffat$xcen
#> [1] 25.5
#>
#> $moffat$ycen
#> [1] 25.5
#>
#> $moffat$mag
#> [1] 0
#>
#> $moffat$fwhm
#> [1] 4.160394
#>
#> $moffat$con
#> [1] 2.460923
#>
#> $moffat$ang
#> [1] 0
#>
#> $moffat$axrat
#> [1] 1
We can now attempt a fit of our target galaxy near the middle of the image. By default this is what is targeted by profitDoFit, but we will specify the loc explicitly to make the process clear. For the first example we will fit just a single Sersic component (the default). The function automatically detects nearby sources and models them a single Sersic profiles to reduce the contamination in our target source.
= profitDoFit(image=image, loc=c(178,178), psf=stars_auto$psf, magzero=30,
example_SS SBdilate=2, Ncomp=1)
#> Running Found2Fit
#> Running ProFound
#> Running Highander
#> Iteration 1
#> CMA 1: -1714.392 29.183 21.93 18.435 1.119 -0.07 45.25 -0.844
#>
#> Laplace's Demon was called on Tue Feb 23 17:52:25 2021
#>
#> Performing initial checks...
#> WARNING: Matrix mask may be rank-deficient.
#> WARNING: Matrix segim may be rank-deficient.
#> WARNING: Matrix psf may be rank-deficient.
#> WARNING: Matrix region may be rank-deficient.
#> WARNING: Matrix calcregion may be rank-deficient.
#> Algorithm: Componentwise Hit-And-Run Metropolis
#>
#> Laplace's Demon is beginning to update...
#> Iteration: 100, Proposal: Componentwise, LP: -1710.3
#>
#> Assessing Stationarity
#> Assessing Thinning and ESS
#> Creating Summaries
#> Creating Output
#>
#> Laplace's Demon has finished.
#> LD Mode 1: -1706.872 29.076 22.053 18.44 1.124 -0.075 44.396 -0.843
#> Iteration 2
#> CMA 2: -1706.721 29.01 22.125 18.441 1.119 -0.067 44.214 -0.841
#>
#> Laplace's Demon was called on Tue Feb 23 17:52:37 2021
#>
#> Performing initial checks...
#> WARNING: Matrix mask may be rank-deficient.
#> WARNING: Matrix segim may be rank-deficient.
#> WARNING: Matrix psf may be rank-deficient.
#> WARNING: Matrix region may be rank-deficient.
#> WARNING: Matrix calcregion may be rank-deficient.
#> Algorithm: Componentwise Hit-And-Run Metropolis
#>
#> Laplace's Demon is beginning to update...
#> Iteration: 100, Proposal: Componentwise, LP: -1709.2
#>
#> Assessing Stationarity
#> Assessing Thinning and ESS
#> Creating Summaries
#> Creating Output
#>
#> Laplace's Demon has finished.
#> LD Mode 2: -1706.496 29.01 22.125 18.441 1.119 -0.067 44.398 -0.841
And check the resulting fit:
profitLikeModel(example_SS$parm, example_SS$Data, makeplots = TRUE, plotchisq = TRUE)
#> Warning in dt(x, tdof): NaNs produced
#> $LP
#> [1] -1706.496
#>
#> $Dev
#> [1] 3412.993
#>
#> $Monitor
#> LL LP tend.elapsed
#> -1706.496 -1706.496 39.296
#>
#> $yhat
#> [1] 1
#>
#> $parm
#> sersic.xcen sersic.ycen sersic.mag sersic.re sersic.nser sersic.ang
#> 29.0100415 22.1249149 18.4406096 1.1188187 -0.0665719 44.3980868
#> sersic.axrat
#> -0.8411036
We can try to improve the fit by doing a 2 component fit. By setting Ncomp=1.5 we will fit the bulge as an unresolved PSF:
= profitDoFit(image=image, loc=c(178,178), psf=stars_auto$psf, magzero=30,
example_Bpsf SBdilate=2, Ncomp=1.5)
#> Running Found2Fit
#> Running ProFound
#> Running Highander
#> Iteration 1
#> CMA 1: -1802.344 28.84 22.351 18.296 1.226 44.447 -0.949 22.922
#>
#> Laplace's Demon was called on Tue Feb 23 17:52:51 2021
#>
#> Performing initial checks...
#> WARNING: Matrix mask may be rank-deficient.
#> WARNING: Matrix segim may be rank-deficient.
#> WARNING: Matrix psf may be rank-deficient.
#> WARNING: Matrix region may be rank-deficient.
#> WARNING: Matrix calcregion may be rank-deficient.
#> Algorithm: Componentwise Hit-And-Run Metropolis
#>
#> Laplace's Demon is beginning to update...
#> Iteration: 100, Proposal: Componentwise, LP: -1793.6
#>
#> Assessing Stationarity
#> Assessing Thinning and ESS
#> Creating Summaries
#> Creating Output
#>
#> Laplace's Demon has finished.
#> LD Mode 1: -1788.907 28.885 22.26 18.283 1.221 44.388 -0.96 25.389
#> Iteration 2
#> CMA 2: -1787.691 28.914 22.246 18.283 1.213 44.299 -0.937 25.91
#>
#> Laplace's Demon was called on Tue Feb 23 17:53:02 2021
#>
#> Performing initial checks...
#> WARNING: Matrix mask may be rank-deficient.
#> WARNING: Matrix segim may be rank-deficient.
#> WARNING: Matrix psf may be rank-deficient.
#> WARNING: Matrix region may be rank-deficient.
#> WARNING: Matrix calcregion may be rank-deficient.
#> Algorithm: Componentwise Hit-And-Run Metropolis
#>
#> Laplace's Demon is beginning to update...
#> Iteration: 100, Proposal: Componentwise, LP: -1788.2
#>
#> Assessing Stationarity
#> Assessing Thinning and ESS
#> Creating Summaries
#> Creating Output
#>
#> Laplace's Demon has finished.
#> LD Mode 2: -1787.525 28.967 22.188 18.285 1.212 44.541 -0.926 27.748
And check the resulting fit:
profitLikeModel(example_Bpsf$parm, example_Bpsf$Data, makeplots = TRUE, plotchisq = TRUE)
#> Warning in dt(x, tdof): NaNs produced
#> $LP
#> [1] -1787.525
#>
#> $Dev
#> [1] 3575.05
#>
#> $Monitor
#> LL LP tend.elapsed
#> -1787.525 -1787.525 64.125
#>
#> $yhat
#> [1] 1
#>
#> $parm
#> sersic.xcen sersic.ycen sersic.mag sersic.re sersic.ang
#> 28.9671096 22.1883413 18.2850905 1.2124503 44.5409697
#> sersic.axrat pointsource.mag
#> -0.9257032 27.7478958
We can try to improve the fit by doing a 2 component fit with a fixed Sersic n=4 bulge.
= profitDoFit(image=image, loc=c(178,178), psf=stars_auto$psf, magzero=30,
example_B4 SBdilate=2, bulge_nser_fit=FALSE, Ncomp=2)
#> Running Found2Fit
#> Running ProFound
#> Running Highander
#> Iteration 1
#> CMA 1: -1714.906 28.952 22.19 20.179 18.483 1.534 1.134 44.61 -0.914
#>
#> Laplace's Demon was called on Tue Feb 23 17:53:20 2021
#>
#> Performing initial checks...
#> WARNING: Matrix mask may be rank-deficient.
#> WARNING: Matrix segim may be rank-deficient.
#> WARNING: Matrix psf may be rank-deficient.
#> WARNING: Matrix region may be rank-deficient.
#> WARNING: Matrix calcregion may be rank-deficient.
#> Algorithm: Componentwise Hit-And-Run Metropolis
#>
#> Laplace's Demon is beginning to update...
#> Iteration: 100, Proposal: Componentwise, LP: -1710.6
#>
#> Assessing Stationarity
#> Assessing Thinning and ESS
#> Creating Summaries
#> Creating Output
#>
#> Laplace's Demon has finished.
#> LD Mode 1: -1709.231 28.995 22.156 21.665 18.456 0.959 1.134 44.362 -0.879
#> Iteration 2
#>
#> Laplace's Demon was called on Tue Feb 23 17:53:42 2021
#>
#> Performing initial checks...
#> WARNING: Matrix mask may be rank-deficient.
#> WARNING: Matrix segim may be rank-deficient.
#> WARNING: Matrix psf may be rank-deficient.
#> WARNING: Matrix region may be rank-deficient.
#> WARNING: Matrix calcregion may be rank-deficient.
#> Algorithm: Componentwise Hit-And-Run Metropolis
#>
#> Laplace's Demon is beginning to update...
#> Iteration: 100, Proposal: Componentwise, LP: -1713.8
#>
#> Assessing Stationarity
#> Assessing Thinning and ESS
#> Creating Summaries
#> Creating Output
#>
#> Laplace's Demon has finished.
And check the resulting fit:
profitLikeModel(example_B4$parm, example_B4$Data, makeplots = TRUE, plotchisq = TRUE)
#> Warning in dt(x, tdof): NaNs produced
#> $LP
#> [1] -1709.231
#>
#> $Dev
#> [1] 3418.462
#>
#> $Monitor
#> LL LP tend.elapsed
#> -1709.231 -1709.231 107.709
#>
#> $yhat
#> [1] 1
#>
#> $parm
#> sersic.xcen1 sersic.ycen1 sersic.mag1 sersic.mag2 sersic.re1
#> 28.9952205 22.1562931 21.6652229 18.4563641 0.9585412
#> sersic.re2 sersic.ang2 sersic.axrat2
#> 1.1338273 44.3618763 -0.8793158