This demo will demonstrate that the different simulation modes for Point Mass lenses are consistent. It may be useful also to look at Pointmass Demo which also includes the closed form simulators in CosmoSim v3.1.
Preparation¶
import pandas as pd
import matplotlib.pyplot as plt
import numpy as np
from PIL import Image
import toml
from copy import deepcopy
import CosmoSim as cs
from CosmoSim import Parameters
from CosmoSim.datagen import SimImage
import CosmoSim.Image as csimg
print( "CosmoSim version:", cs.__version__ )CosmoSim version: 3.2.2b2
Raytrace¶
We can define the configuration as a dict using the nested (TOML) structure.
cfg = { 'simulator' : { "model" : "Raytrace", "nterms" : 10, "cropsize" : 256 }
, 'lens': {
'mode' : "PM",
'einsteinradius': 46 }
, 'source': {
'mode': 'SersicSphere',
'sigma': 20,
'theta': 45,
'luminosity' : 20,
'position': 'cartesian'}
, 'position': {'x': 11.01, 'y': 0.31}
}
param = Parameters(cfg)raysim = SimImage( param, verbose=0 )
rayim = raysim.getImage()
csimg.imshow( rayim, title="Raytrace Point Mass")
Roulette¶
param["simulator"]["model"] = "Roulette"
rousim = SimImage( param, verbose=2 )
rouim = rousim.getImage()
csimg.imageCompare( rouim, rayim, "Roulette model", 'Raytrace')[GenericSim] init (verbose=2) ...
[SimImage] init (verbose=2) ...
[getSource] src=SersicSphere, ltprf0=None, verbose=2
[getSource] Lightprofile: SersicSphere LightProfileSpec.Sersic
[getSource] mode=SourceSpec.Sphere, ltprf=LightProfileSpec.Sersic
[getSource] Spherical Source - n_sersic=4.0, luminosity=20.0
[SphericalSource] constructor done
getSource() returns
[initSim] XY 11.01 0.31
[Source] Constructor
[SphericalSource] SERSIC
setFile /home/runner/.local/lib/python3.14/site-packages/CosmoSim/pm50.txt
[PsiFunctionLens.initAlphasBetas] Amplitudes file /home/runner/.local/lib/python3.14/site-packages/CosmoSim/pm50.txt
[initAlphasBetas] opened file /home/runner/.local/lib/python3.14/site-packages/CosmoSim/pm50.txt
Instantiating RouletteModel ...
[SimulatorModel::setSource] setting source
[SimulatorModel::setNterms] 10 -> 10
[SimulatorModel::update] Lens: PointMass
[PointMass::getXi] [11.01, 0.31] -> [51.8151, 1.45892]
[setNu] etaOffset set to zero.
[SimulatorModel::update] Done updateApparentAbs()
[SimulatorModel::update] thread section
[Source::getImage()]
[SimulatorModel::updateInner()] eta=[11.01, 0.31]
[SimulatorModel::updateInner()] xi=[51.8151, 1.45892]; eta=[11.01, 0.31]; etaOffset=[0, 0]
[SimulatorModel::updateInner()] nu=[51.8151, 1.45892]
[calculateAlphaBeta] [[51.8151, 1.45892]] ...
[PsiFunctionLens.calculateAlphaBeta()] 10; 46 - [51.8151, 1.45892]
[SimulatorModel::parallelDistort] 4 threads (maskMode=0)
[SimulatorModel::parallelDistort] etaOffset=[0, 0]; nu=[51.8151, 1.45892]; eta[11.01, 0.31]; referenceXi=[51.8151, 1.45892]
[SimulatorModel] No mask
[SimulatorModel::parallelDistort] lower=0; rng=512; rng1=128
[SimulatorModel::distort] refXi=[51.8151, 1.45892] begin=0; end=128
[SimulatorModel::distort] refXi=[51.8151, 1.45892] begin=384; end=512
[SimulatorModel::distort] refXi=[51.8151, 1.45892] begin=128; end=256
[SimulatorModel::distort] refXi=[51.8151, 1.45892] begin=256; end=384
[getDistortedPos(2)] nu=[-1.5002, -0.31] theta=-2.08032 r=1.67119 (PointMass) [1.67119, -2.08032]->[-1.5002, -0.31]
[getDistortedPos(2)] nu=[0.297692, -0.31] theta=-1.44476 r=1.47058 (PointMass) [1.47058, -1.44476]->[0.297692, -0.31]
[getDistortedPos(2)] nu=[2.06547, -0.31] theta=-0.88869 r=1.87945 (PointMass) [1.87945, -0.88869]->[2.06547, -0.31]
[getDistortedPos(2)] nu=[-1.35713, 0.257816] theta=2.05733 r=1.74338 (PointMass) [1.74338, 2.05733]->[-1.35713, 0.257816]
[getDistortedPos(2)] nu=[0.432683, 0.350155] theta=1.45141 r=1.55213 (PointMass) [1.55213, 1.45141]->[0.432683, 0.350155]
[getDistortedPos(2)] nu=[2.19298, 0.437339] theta=0.915339 r=1.94392 (PointMass) [1.94392, 0.915339]->[2.19298, 0.437339]
[getDistortedPos(2)] nu=[-3.2624, -9.58466e-05] theta=2.85189 r=1.89407 (PointMass) [1.89407, 2.85189]->[-3.2624, -9.58466e-05]
[getDistortedPos(2)] nu=[-1.43649, 0.0654319] theta=2.55558 r=0.978375 (PointMass) [0.978375, 2.55558]->[-1.43649, 0.0654319]
[getDistortedPos(2)] nu=[0.357799, 0.127223] theta=1.24158 r=0.57179 (PointMass) [0.57179, 1.24158]->[0.357799, 0.127223]
[getDistortedPos(2)] nu=[2.12224, 0.185556] theta=0.42838 r=1.30256 (PointMass) [1.30256, 0.42838]->[2.12224, 0.185556]
[getDistortedPos(2)] nu=[-3.31308, -0.156062] theta=-2.89395 r=1.87225 (PointMass) [1.87225, -2.89395]->[-3.31308, -0.156062]
[getDistortedPos(2)] nu=[-1.48425, -0.123221] theta=-2.62883 r=0.935444 (PointMass) [0.935444, -2.62883]->[-1.48425, -0.123221]
[getDistortedPos(2)] nu=[0.312736, -0.0922551] theta=-1.18786 r=0.494753 (PointMass) [0.494753, -1.18786]->[0.312736, -0.0922551]
[getDistortedPos(2)] nu=[2.07968, -0.0630249] theta=-0.369526 r=1.27063 (PointMass) [1.27063, -0.369526]->[2.07968, -0.0630249]
Time to update(): 196 milliseconds
[Simulator] Centre Point (48.08,0.12) (Centre of Luminence in Planar Co-ordinates)
[getImage] centring from parameters: None
[crop] cropsize=256
This looks perfect inside the convergence ring, as it should.
Sampling¶
param["simulator"]["sampled"] = True
imsim05 = SimImage( param, verbose=0 )
im05 = imsim05.getImage()
csimg.imageCompare( rayim, im05, "Raytrace", 'Sampled Roulette')
csimg.imageCompare( rouim, im05, "Roulette", 'Sampled Roulette')[SimulatorModel::getDistorted()]
[setDebug] 0
