[1]:
import numpy as np
import clusterps
Setup
[2]:
kT = 8.54
boxsize = 60. # arcmin
nn=200
z_coma = 0.023
[3]:
cps = clusterps.ClusterPS(redshift=z_coma, kT=kT, boxsize=boxsize, npix=nn)
No cosmology provided, will default to Planck15
Luminosity distance to the source: 103.573 Mpc
At the redshift of the source 1 arcmin is 28.7887 kpc
We will use a box of 200 pixels on a side, corresponding to a physical length of 1727.32 kpc and a resolution of 0.3 kpc
Given the provided temperature of 8.54 keV, we will compute the normalization of the power spectrum in terms of the Mach number respective to the sound speed cs=1495.01 km/s in the medium
[4]:
rc = 11.7 * cps.cosmo.kpc_proper_per_arcmin(z_coma).value
beta = 0.765
emprof = lambda x: (1 + (x/rc)**2) ** (-3. * beta)
[5]:
cps.EMmodel(emprof)
Projected EM model saved to file PS/emmodel.fits
[6]:
cps.SetRegion('region_def.reg')
box(92.344781
91.344781 104.35994 5.0998889 16.000007
2
box(97.252752
96.252752 105.71213 5.0998333 16.000007
1
box(93.696987
92.696987 99.45196 5.0998333 16.000007
3
box(98.604947
97.604947 100.85424 5.0998333 16.000007
4
box(95.520003
94.520003 79.189341 5.0998889 44.000007
8
box(92.00399
91.00399 82.875241 5.0998333 44.000007
6
box(99.226012
98.226012 82.745121 5.0998333 44.000007
7
box(95.700183
94.700183 86.41088 5.0998333 44.000007
5
Region definition map saved to file PS/regions.fits
[7]:
cps.SetRegion('region_def_vdisp.reg', outfile='PS/regions_vdisp.fits', observable='disp')
box(92.344781
91.344781 104.35994 5.0998889 16.000007
2
box(97.252752
96.252752 105.71213 5.0998333 16.000007
1
box(93.696987
92.696987 99.45196 5.0998333 16.000007
3
box(98.604947
97.604947 100.85424 5.0998333 16.000007
4
box(95.640336
94.640336 82.804208 10.232681 44.000007
5
Region definition map saved to file PS/regions_vdisp.fits
[8]:
velocities = np.array([-415., -400., -410., -570., -770., -690., -800., -690.])
vb_errors = np.array([30., 30., 30., 30., 60., 60., 60., 60.])
vdisp = np.array([200., 250., 145., 210., 200.])
vdisp_errors = np.array([25., 25., 25., 25., 23.])
[9]:
cps.SetBulkVelocities(velocities, vb_errors)
cps.SetVelocityDispersions(vdisp, vdisp_errors)
[10]:
cps.PSF(34.)
cps.psf
[10]:
1.8888888888888888
[11]:
def fturb(k, pars):
slope = pars[0]
kinj = 10 ** pars[1]
return (1 + (k/kinj)**2 ) ** (-slope / 2.)
mod = clusterps.Model(model=fturb, npar=3)
[13]:
mod.SetUniformPrior(0, 0., 2.)
mod.SetUniformPrior(1, 2., 8.)
mod.SetUniformPrior(2, -4, -2)
[14]:
mod.priorlist
[14]:
[BoxUniform(Uniform(low: tensor([0.]), high: tensor([2.])), 1),
BoxUniform(Uniform(low: tensor([2.]), high: tensor([8.])), 1),
BoxUniform(Uniform(low: tensor([-4.]), high: tensor([-2.])), 1)]
[15]:
mod.SetPrior()
[16]:
mod.prior
[16]:
MultipleIndependent()
Run
[21]:
clusterps.sbi.SBISetup(cps, mod)
Initializing 3D k-space ...
Initializing fluctuations in the k-space ...
[24]:
clusterps.sbi.SimulateForSBI(cps, mod)
Neural network successfully converged after 202 epochs.
[25]:
clusterps.sbi.Infer(cps, mod)
[26]:
mod.chains
[26]:
array([[ 0.7164708 , 7.6973867 , -3.3626423 ],
[ 0.7531437 , 6.697179 , -3.4504263 ],
[ 0.84699446, 4.3825374 , -3.5711029 ],
...,
[ 0.60852665, 3.8092232 , -3.919593 ],
[ 1.0892688 , 6.844991 , -3.1631675 ],
[ 0.8086788 , 7.0572896 , -3.0737991 ]], dtype=float32)
[29]:
import pandas as pd
from chainconsumer import Chain, ChainConsumer
c = ChainConsumer()
np_samples = mod.chains
df = pd.DataFrame({
'$\mathcal{M}_{3D}$': np_samples[:,0],
'$\\alpha$': np_samples[:,1],
'$k_{inj}$': np_samples[:,2],
})
c.add_chain(Chain(samples=df, name='Coma 2 pointings'))
fig = c.plotter.plot()
Parameter $\alpha$ in chain Coma 2 pointings is not constrained
Reload
[30]:
cps.Reload(infile='PS/sbi_model.pkl', Model=mod, chains='PS/sbi_inference.npy')
Initializing 3D k-space ...
Initializing fluctuations in the k-space ...
/Users/deckert/opt/anaconda3/envs/sbi_env/lib/python3.11/site-packages/torch/storage.py:414: FutureWarning: You are using `torch.load` with `weights_only=False` (the current default value), which uses the default pickle module implicitly. It is possible to construct malicious pickle data which will execute arbitrary code during unpickling (See https://github.com/pytorch/pytorch/blob/main/SECURITY.md#untrusted-models for more details). In a future release, the default value for `weights_only` will be flipped to `True`. This limits the functions that could be executed during unpickling. Arbitrary objects will no longer be allowed to be loaded via this mode unless they are explicitly allowlisted by the user via `torch.serialization.add_safe_globals`. We recommend you start setting `weights_only=True` for any use case where you don't have full control of the loaded file. Please open an issue on GitHub for any issues related to this experimental feature.
return torch.load(io.BytesIO(b))
[37]:
import torch
cps.simulator(torch.tensor([0.5, 11./3., -3.5]))
[37]:
tensor([327.6991, 256.5038, 317.1894, 409.2164, 354.0026, 186.1683, 193.2260,
64.4723, 167.6507, 172.8217, 103.1162, 132.5215, 234.1563],
dtype=torch.float64)
[31]:
clusterps.sbi.Infer(cps, mod)
TARP tests
[32]:
import sbi
sbi.__version__
[32]:
'0.23.3'
[29]:
thetas = mod.prior.sample((1000,))
xs = cps.simulator(thetas)
[43]:
from tqdm import tqdm
import torch
nsim = 1000
xs = torch.zeros(nsim, cps.nval)
for i in tqdm(range(nsim)):
xs[i] = cps.simulator(thetas[i])
100%|███████████████████████████████████████| 1000/1000 [10:03<00:00, 1.66it/s]
[47]:
outpar = np.empty((1000, mod.npar))
nsamples = 10000
for i in tqdm(range(1000)):
samples = cps.inferencemodel.sample((nsamples,), x=xs[i])
outpar[i] = np.median(samples, axis=0)
0%| | 0/1000 [00:00<?, ?it/s]
0%| | 1/1000 [00:00<02:07, 7.85it/s]
0%|▏ | 3/1000 [00:00<01:36, 10.33it/s]
0%|▏ | 5/1000 [00:00<02:09, 7.68it/s]
1%|▎ | 6/1000 [00:00<02:08, 7.71it/s]
1%|▎ | 8/1000 [00:00<01:59, 8.32it/s]
1%|▍ | 10/1000 [00:01<01:46, 9.31it/s]
1%|▍ | 11/1000 [00:01<01:44, 9.44it/s]
1%|▌ | 13/1000 [00:01<01:42, 9.65it/s]
1%|▌ | 14/1000 [00:01<01:50, 8.90it/s]
2%|▋ | 16/1000 [00:01<01:43, 9.52it/s]
2%|▋ | 18/1000 [00:01<01:37, 10.05it/s]
2%|▊ | 19/1000 [00:02<01:40, 9.80it/s]
2%|▊ | 21/1000 [00:02<01:36, 10.13it/s]
2%|▉ | 23/1000 [00:02<01:35, 10.26it/s]
2%|█ | 25/1000 [00:02<01:30, 10.82it/s]
3%|█ | 27/1000 [00:02<01:28, 10.99it/s]
3%|█▏ | 29/1000 [00:02<01:26, 11.18it/s]
3%|█▎ | 31/1000 [00:03<01:30, 10.72it/s]
3%|█▎ | 33/1000 [00:03<01:30, 10.72it/s]
4%|█▍ | 35/1000 [00:03<01:31, 10.55it/s]
4%|█▌ | 37/1000 [00:03<01:27, 11.07it/s]
4%|█▌ | 39/1000 [00:03<01:23, 11.55it/s]
4%|█▋ | 41/1000 [00:04<01:22, 11.61it/s]
4%|█▊ | 43/1000 [00:04<01:24, 11.36it/s]
4%|█▊ | 45/1000 [00:04<01:25, 11.23it/s]
5%|█▉ | 47/1000 [00:04<01:22, 11.60it/s]
5%|██ | 49/1000 [00:04<01:21, 11.71it/s]
5%|██ | 51/1000 [00:04<01:34, 10.06it/s]
5%|██▏ | 53/1000 [00:05<01:34, 10.02it/s]
6%|██▎ | 55/1000 [00:05<01:41, 9.30it/s]
6%|██▎ | 57/1000 [00:05<01:35, 9.89it/s]
6%|██▍ | 59/1000 [00:05<01:28, 10.62it/s]
6%|██▌ | 61/1000 [00:05<01:31, 10.24it/s]
6%|██▌ | 63/1000 [00:06<01:28, 10.63it/s]
6%|██▋ | 65/1000 [00:06<01:44, 8.98it/s]
7%|██▋ | 67/1000 [00:06<01:45, 8.82it/s]
7%|██▊ | 69/1000 [00:06<01:42, 9.06it/s]
7%|██▊ | 70/1000 [00:07<01:44, 8.88it/s]
7%|██▉ | 72/1000 [00:07<01:42, 9.08it/s]
7%|███ | 74/1000 [00:07<01:35, 9.73it/s]
8%|███ | 76/1000 [00:07<01:31, 10.05it/s]
8%|███▏ | 78/1000 [00:07<01:30, 10.15it/s]
8%|███▎ | 80/1000 [00:08<01:43, 8.85it/s]
8%|███▎ | 81/1000 [00:08<01:50, 8.30it/s]
8%|███▎ | 82/1000 [00:08<01:54, 8.00it/s]
8%|███▍ | 84/1000 [00:08<01:42, 8.96it/s]
9%|███▌ | 86/1000 [00:08<01:31, 9.96it/s]
9%|███▌ | 88/1000 [00:08<01:25, 10.66it/s]
9%|███▋ | 90/1000 [00:09<01:29, 10.14it/s]
9%|███▊ | 92/1000 [00:09<01:44, 8.70it/s]
9%|███▊ | 94/1000 [00:09<01:38, 9.20it/s]
10%|███▉ | 95/1000 [00:09<01:43, 8.71it/s]
10%|███▉ | 97/1000 [00:09<01:35, 9.45it/s]
10%|████ | 99/1000 [00:10<01:31, 9.90it/s]
10%|████ | 101/1000 [00:10<01:25, 10.48it/s]
10%|████ | 103/1000 [00:10<01:21, 11.05it/s]
10%|████▏ | 105/1000 [00:10<01:18, 11.46it/s]
11%|████▎ | 107/1000 [00:10<01:17, 11.54it/s]
11%|████▎ | 109/1000 [00:10<01:22, 10.79it/s]
11%|████▍ | 111/1000 [00:11<01:19, 11.18it/s]
11%|████▌ | 113/1000 [00:11<01:18, 11.35it/s]
12%|████▌ | 115/1000 [00:11<01:19, 11.09it/s]
12%|████▋ | 117/1000 [00:11<01:23, 10.61it/s]
12%|████▊ | 119/1000 [00:11<01:21, 10.86it/s]
12%|████▊ | 121/1000 [00:12<01:17, 11.28it/s]
12%|████▉ | 123/1000 [00:12<01:22, 10.64it/s]
12%|█████ | 125/1000 [00:12<01:19, 11.02it/s]
13%|█████ | 127/1000 [00:12<01:19, 10.98it/s]
13%|█████▏ | 129/1000 [00:12<01:20, 10.87it/s]
13%|█████▏ | 131/1000 [00:12<01:18, 11.13it/s]
13%|█████▎ | 133/1000 [00:13<01:18, 10.98it/s]
14%|█████▍ | 135/1000 [00:13<01:26, 10.05it/s]
14%|█████▍ | 137/1000 [00:13<01:29, 9.66it/s]
14%|█████▌ | 139/1000 [00:13<01:24, 10.19it/s]
14%|█████▋ | 141/1000 [00:13<01:21, 10.50it/s]
14%|█████▋ | 143/1000 [00:14<01:18, 10.96it/s]
14%|█████▊ | 145/1000 [00:14<01:15, 11.39it/s]
15%|█████▉ | 147/1000 [00:14<01:18, 10.91it/s]
15%|█████▉ | 149/1000 [00:14<01:17, 10.97it/s]
15%|██████ | 151/1000 [00:14<01:12, 11.66it/s]
15%|██████ | 153/1000 [00:14<01:12, 11.67it/s]
16%|██████▏ | 155/1000 [00:15<01:13, 11.42it/s]
16%|██████▎ | 157/1000 [00:15<01:15, 11.21it/s]
16%|██████▎ | 159/1000 [00:15<01:12, 11.56it/s]
16%|██████▍ | 161/1000 [00:15<01:17, 10.79it/s]
16%|██████▌ | 163/1000 [00:15<01:21, 10.25it/s]
16%|██████▌ | 165/1000 [00:16<01:21, 10.30it/s]
17%|██████▋ | 167/1000 [00:16<01:16, 10.86it/s]
17%|██████▊ | 169/1000 [00:16<01:25, 9.75it/s]
17%|██████▊ | 171/1000 [00:17<02:02, 6.78it/s]
17%|██████▉ | 172/1000 [00:17<01:56, 7.08it/s]
17%|██████▉ | 173/1000 [00:17<01:52, 7.32it/s]
17%|██████▉ | 174/1000 [00:17<01:48, 7.58it/s]
18%|███████ | 176/1000 [00:17<01:40, 8.19it/s]
18%|███████ | 178/1000 [00:17<01:32, 8.84it/s]
18%|███████▏ | 179/1000 [00:17<01:33, 8.75it/s]
18%|███████▏ | 181/1000 [00:18<01:31, 8.98it/s]
18%|███████▎ | 183/1000 [00:18<01:24, 9.66it/s]
18%|███████▍ | 185/1000 [00:18<01:21, 10.01it/s]
19%|███████▍ | 187/1000 [00:18<01:18, 10.31it/s]
19%|███████▌ | 189/1000 [00:18<01:14, 10.94it/s]
19%|███████▋ | 191/1000 [00:19<01:22, 9.84it/s]
19%|███████▋ | 193/1000 [00:19<01:40, 8.03it/s]
19%|███████▊ | 194/1000 [00:19<01:36, 8.32it/s]
20%|███████▊ | 195/1000 [00:19<01:34, 8.53it/s]
20%|███████▊ | 196/1000 [00:19<01:33, 8.63it/s]
20%|███████▉ | 198/1000 [00:19<01:26, 9.28it/s]
20%|████████ | 200/1000 [00:20<01:18, 10.19it/s]
20%|████████ | 202/1000 [00:20<01:24, 9.41it/s]
20%|████████ | 203/1000 [00:20<01:26, 9.17it/s]
20%|████████▏ | 204/1000 [00:20<01:27, 9.10it/s]
21%|████████▏ | 206/1000 [00:20<01:22, 9.64it/s]
21%|████████▎ | 207/1000 [00:20<01:22, 9.63it/s]
21%|████████▎ | 208/1000 [00:20<01:23, 9.46it/s]
21%|████████▍ | 210/1000 [00:21<01:18, 10.07it/s]
21%|████████▍ | 211/1000 [00:21<01:19, 9.98it/s]
21%|████████▍ | 212/1000 [00:21<01:24, 9.35it/s]
21%|████████▌ | 213/1000 [00:21<01:28, 8.94it/s]
21%|████████▌ | 214/1000 [00:21<01:27, 8.93it/s]
22%|████████▌ | 215/1000 [00:21<01:30, 8.71it/s]
22%|████████▋ | 216/1000 [00:21<01:36, 8.12it/s]
22%|████████▋ | 217/1000 [00:22<01:40, 7.81it/s]
22%|████████▋ | 218/1000 [00:22<01:39, 7.88it/s]
22%|████████▊ | 219/1000 [00:22<01:40, 7.80it/s]
22%|████████▊ | 221/1000 [00:22<01:24, 9.23it/s]
22%|████████▉ | 223/1000 [00:22<01:16, 10.19it/s]
22%|█████████ | 225/1000 [00:22<01:10, 11.00it/s]
23%|█████████ | 227/1000 [00:23<01:20, 9.56it/s]
23%|█████████ | 228/1000 [00:23<01:24, 9.10it/s]
23%|█████████▏ | 230/1000 [00:23<01:17, 9.97it/s]
23%|█████████▎ | 232/1000 [00:23<01:13, 10.46it/s]
23%|█████████▎ | 234/1000 [00:23<01:12, 10.63it/s]
24%|█████████▍ | 236/1000 [00:23<01:11, 10.75it/s]
24%|█████████▌ | 238/1000 [00:24<01:17, 9.77it/s]
24%|█████████▌ | 239/1000 [00:24<01:18, 9.72it/s]
24%|█████████▋ | 241/1000 [00:24<01:13, 10.31it/s]
24%|█████████▋ | 243/1000 [00:24<01:12, 10.48it/s]
24%|█████████▊ | 245/1000 [00:24<01:07, 11.21it/s]
25%|█████████▉ | 247/1000 [00:24<01:05, 11.46it/s]
25%|█████████▉ | 249/1000 [00:25<01:04, 11.61it/s]
25%|██████████ | 251/1000 [00:25<01:03, 11.84it/s]
25%|██████████ | 253/1000 [00:25<01:06, 11.27it/s]
26%|██████████▏ | 255/1000 [00:25<01:05, 11.35it/s]
26%|██████████▎ | 257/1000 [00:25<01:03, 11.62it/s]
26%|██████████▎ | 259/1000 [00:25<01:03, 11.59it/s]
26%|██████████▍ | 261/1000 [00:26<01:05, 11.25it/s]
26%|██████████▌ | 263/1000 [00:26<01:21, 9.06it/s]
26%|██████████▌ | 264/1000 [00:26<01:26, 8.55it/s]
27%|██████████▋ | 266/1000 [00:26<01:19, 9.21it/s]
27%|██████████▋ | 268/1000 [00:26<01:14, 9.86it/s]
27%|██████████▊ | 270/1000 [00:27<01:11, 10.16it/s]
27%|██████████▉ | 272/1000 [00:27<01:12, 10.00it/s]
27%|██████████▉ | 274/1000 [00:27<01:12, 10.02it/s]
28%|███████████ | 276/1000 [00:27<01:11, 10.17it/s]
28%|███████████ | 278/1000 [00:27<01:10, 10.24it/s]
28%|███████████▏ | 280/1000 [00:28<01:07, 10.74it/s]
28%|███████████▎ | 282/1000 [00:28<01:02, 11.40it/s]
28%|███████████▎ | 284/1000 [00:28<01:01, 11.69it/s]
29%|███████████▍ | 286/1000 [00:28<01:03, 11.17it/s]
29%|███████████▌ | 288/1000 [00:28<01:01, 11.59it/s]
29%|███████████▌ | 290/1000 [00:28<01:00, 11.69it/s]
29%|███████████▋ | 292/1000 [00:29<01:00, 11.70it/s]
29%|███████████▊ | 294/1000 [00:29<01:00, 11.62it/s]
30%|███████████▊ | 296/1000 [00:29<00:58, 12.11it/s]
30%|███████████▉ | 298/1000 [00:29<01:01, 11.40it/s]
30%|████████████ | 300/1000 [00:29<01:07, 10.36it/s]
30%|████████████ | 302/1000 [00:30<01:07, 10.33it/s]
30%|████████████▏ | 304/1000 [00:30<01:05, 10.67it/s]
31%|████████████▏ | 306/1000 [00:30<01:18, 8.88it/s]
31%|████████████▎ | 308/1000 [00:30<01:12, 9.54it/s]
31%|████████████▍ | 310/1000 [00:30<01:10, 9.85it/s]
31%|████████████▍ | 312/1000 [00:31<01:09, 9.94it/s]
31%|████████████▌ | 314/1000 [00:31<01:08, 10.08it/s]
32%|████████████▋ | 316/1000 [00:31<01:05, 10.37it/s]
32%|████████████▋ | 318/1000 [00:31<01:04, 10.54it/s]
32%|████████████▊ | 320/1000 [00:31<01:03, 10.79it/s]
32%|████████████▉ | 322/1000 [00:31<00:59, 11.39it/s]
32%|████████████▉ | 324/1000 [00:32<00:58, 11.60it/s]
33%|█████████████ | 326/1000 [00:32<00:57, 11.75it/s]
33%|█████████████ | 328/1000 [00:32<00:58, 11.45it/s]
33%|█████████████▏ | 330/1000 [00:32<00:58, 11.44it/s]
33%|█████████████▎ | 332/1000 [00:32<00:56, 11.90it/s]
33%|█████████████▎ | 334/1000 [00:32<00:54, 12.15it/s]
34%|█████████████▍ | 336/1000 [00:33<00:53, 12.43it/s]
34%|█████████████▌ | 338/1000 [00:33<00:53, 12.46it/s]
34%|█████████████▌ | 340/1000 [00:33<00:54, 12.18it/s]
34%|█████████████▋ | 342/1000 [00:33<00:58, 11.33it/s]
34%|█████████████▊ | 344/1000 [00:33<00:58, 11.18it/s]
35%|█████████████▊ | 346/1000 [00:33<00:56, 11.53it/s]
35%|█████████████▉ | 348/1000 [00:34<00:54, 11.87it/s]
35%|██████████████ | 350/1000 [00:34<00:54, 11.85it/s]
35%|██████████████ | 352/1000 [00:34<00:53, 12.16it/s]
35%|██████████████▏ | 354/1000 [00:34<00:53, 12.05it/s]
36%|██████████████▏ | 356/1000 [00:34<00:52, 12.26it/s]
36%|██████████████▎ | 358/1000 [00:34<00:53, 11.91it/s]
36%|██████████████▍ | 360/1000 [00:35<01:08, 9.39it/s]
36%|██████████████▍ | 362/1000 [00:35<01:11, 8.94it/s]
36%|██████████████▌ | 363/1000 [00:35<01:14, 8.52it/s]
36%|██████████████▌ | 364/1000 [00:35<01:13, 8.63it/s]
36%|██████████████▌ | 365/1000 [00:35<01:18, 8.11it/s]
37%|██████████████▋ | 366/1000 [00:36<01:14, 8.50it/s]
37%|██████████████▋ | 368/1000 [00:36<01:10, 9.03it/s]
37%|██████████████▊ | 370/1000 [00:36<01:02, 10.15it/s]
37%|██████████████▉ | 372/1000 [00:36<01:00, 10.32it/s]
37%|██████████████▉ | 374/1000 [00:36<00:59, 10.50it/s]
38%|███████████████ | 376/1000 [00:37<01:05, 9.57it/s]
38%|███████████████ | 378/1000 [00:37<01:05, 9.53it/s]
38%|███████████████▏ | 379/1000 [00:37<01:05, 9.51it/s]
38%|███████████████▏ | 380/1000 [00:37<01:04, 9.61it/s]
38%|███████████████▏ | 381/1000 [00:37<01:04, 9.55it/s]
38%|███████████████▎ | 382/1000 [00:37<01:04, 9.58it/s]
38%|███████████████▎ | 384/1000 [00:37<01:02, 9.84it/s]
38%|███████████████▍ | 385/1000 [00:37<01:09, 8.86it/s]
39%|███████████████▍ | 386/1000 [00:38<01:11, 8.59it/s]
39%|███████████████▌ | 388/1000 [00:38<01:05, 9.41it/s]
39%|███████████████▌ | 390/1000 [00:38<00:59, 10.29it/s]
39%|███████████████▋ | 392/1000 [00:38<01:04, 9.49it/s]
39%|███████████████▊ | 394/1000 [00:38<01:07, 9.00it/s]
40%|███████████████▊ | 395/1000 [00:39<01:08, 8.89it/s]
40%|███████████████▊ | 396/1000 [00:39<01:19, 7.63it/s]
40%|███████████████▉ | 398/1000 [00:39<01:09, 8.69it/s]
40%|████████████████ | 400/1000 [00:39<01:03, 9.52it/s]
40%|████████████████ | 402/1000 [00:39<01:01, 9.74it/s]
40%|████████████████▏ | 404/1000 [00:40<01:02, 9.61it/s]
41%|████████████████▏ | 406/1000 [00:40<00:59, 9.93it/s]
41%|████████████████▎ | 408/1000 [00:40<01:01, 9.69it/s]
41%|████████████████▍ | 410/1000 [00:40<00:59, 9.86it/s]
41%|████████████████▍ | 411/1000 [00:40<01:00, 9.67it/s]
41%|████████████████▍ | 412/1000 [00:40<01:01, 9.55it/s]
41%|████████████████▌ | 414/1000 [00:41<00:56, 10.31it/s]
42%|████████████████▋ | 416/1000 [00:41<00:52, 11.20it/s]
42%|████████████████▋ | 418/1000 [00:41<00:51, 11.28it/s]
42%|████████████████▊ | 420/1000 [00:41<00:51, 11.33it/s]
42%|████████████████▉ | 422/1000 [00:41<00:50, 11.54it/s]
42%|████████████████▉ | 424/1000 [00:41<00:48, 11.77it/s]
43%|█████████████████ | 426/1000 [00:42<00:52, 11.02it/s]
43%|█████████████████ | 428/1000 [00:42<00:52, 10.87it/s]
43%|█████████████████▏ | 430/1000 [00:42<00:51, 11.10it/s]
43%|█████████████████▎ | 432/1000 [00:42<00:49, 11.46it/s]
43%|█████████████████▎ | 434/1000 [00:42<00:49, 11.46it/s]
44%|█████████████████▍ | 436/1000 [00:43<00:54, 10.28it/s]
44%|█████████████████▌ | 438/1000 [00:43<00:56, 9.90it/s]
44%|█████████████████▌ | 440/1000 [00:43<01:00, 9.26it/s]
44%|█████████████████▋ | 441/1000 [00:43<01:01, 9.15it/s]
44%|█████████████████▋ | 442/1000 [00:43<01:01, 9.01it/s]
44%|█████████████████▋ | 443/1000 [00:43<01:14, 7.49it/s]
44%|█████████████████▊ | 444/1000 [00:44<01:19, 6.98it/s]
44%|█████████████████▊ | 445/1000 [00:44<01:53, 4.88it/s]
45%|█████████████████▊ | 446/1000 [00:44<01:47, 5.17it/s]
45%|█████████████████▉ | 447/1000 [00:44<01:46, 5.18it/s]
45%|█████████████████▉ | 448/1000 [00:45<01:46, 5.17it/s]
45%|█████████████████▉ | 449/1000 [00:45<01:38, 5.61it/s]
45%|██████████████████ | 450/1000 [00:45<01:30, 6.05it/s]
45%|██████████████████ | 451/1000 [00:45<01:27, 6.31it/s]
45%|██████████████████ | 453/1000 [00:45<01:06, 8.19it/s]
46%|██████████████████▏ | 455/1000 [00:45<00:59, 9.09it/s]
46%|██████████████████▎ | 457/1000 [00:45<00:52, 10.30it/s]
46%|██████████████████▎ | 459/1000 [00:46<00:49, 10.90it/s]
46%|██████████████████▍ | 461/1000 [00:46<00:47, 11.37it/s]
46%|██████████████████▌ | 463/1000 [00:46<00:45, 11.88it/s]
46%|██████████████████▌ | 465/1000 [00:46<00:42, 12.63it/s]
47%|██████████████████▋ | 467/1000 [00:46<00:42, 12.55it/s]
47%|██████████████████▊ | 469/1000 [00:46<00:46, 11.31it/s]
47%|██████████████████▊ | 471/1000 [00:47<00:45, 11.58it/s]
47%|██████████████████▉ | 473/1000 [00:47<00:46, 11.39it/s]
48%|███████████████████ | 475/1000 [00:47<00:44, 11.77it/s]
48%|███████████████████ | 477/1000 [00:47<00:41, 12.47it/s]
48%|███████████████████▏ | 479/1000 [00:47<00:44, 11.73it/s]
48%|███████████████████▏ | 481/1000 [00:47<00:47, 10.87it/s]
48%|███████████████████▎ | 483/1000 [00:48<00:47, 10.78it/s]
48%|███████████████████▍ | 485/1000 [00:48<00:45, 11.20it/s]
49%|███████████████████▍ | 487/1000 [00:48<00:44, 11.65it/s]
49%|███████████████████▌ | 489/1000 [00:48<00:41, 12.29it/s]
49%|███████████████████▋ | 491/1000 [00:48<00:40, 12.70it/s]
49%|███████████████████▋ | 493/1000 [00:48<00:40, 12.42it/s]
50%|███████████████████▊ | 495/1000 [00:49<00:40, 12.56it/s]
50%|███████████████████▉ | 497/1000 [00:49<00:39, 12.81it/s]
50%|███████████████████▉ | 499/1000 [00:49<00:39, 12.64it/s]
50%|████████████████████ | 501/1000 [00:49<00:50, 9.83it/s]
50%|████████████████████ | 503/1000 [00:49<00:51, 9.68it/s]
50%|████████████████████▏ | 505/1000 [00:50<00:51, 9.53it/s]
51%|████████████████████▎ | 507/1000 [00:50<00:50, 9.81it/s]
51%|████████████████████▎ | 509/1000 [00:50<00:50, 9.82it/s]
51%|████████████████████▍ | 511/1000 [00:50<00:47, 10.36it/s]
51%|████████████████████▌ | 513/1000 [00:50<00:50, 9.74it/s]
52%|████████████████████▌ | 515/1000 [00:51<00:48, 10.03it/s]
52%|████████████████████▋ | 517/1000 [00:51<00:45, 10.54it/s]
52%|████████████████████▊ | 519/1000 [00:51<00:43, 11.08it/s]
52%|████████████████████▊ | 521/1000 [00:51<00:40, 11.73it/s]
52%|████████████████████▉ | 523/1000 [00:51<00:40, 11.73it/s]
52%|█████████████████████ | 525/1000 [00:51<00:39, 12.13it/s]
53%|█████████████████████ | 527/1000 [00:52<00:40, 11.54it/s]
53%|█████████████████████▏ | 529/1000 [00:52<00:42, 11.06it/s]
53%|█████████████████████▏ | 531/1000 [00:52<00:43, 10.70it/s]
53%|█████████████████████▎ | 533/1000 [00:52<00:41, 11.21it/s]
54%|█████████████████████▍ | 535/1000 [00:52<00:38, 11.97it/s]
54%|█████████████████████▍ | 537/1000 [00:52<00:40, 11.40it/s]
54%|█████████████████████▌ | 539/1000 [00:53<00:38, 11.82it/s]
54%|█████████████████████▋ | 541/1000 [00:53<00:38, 11.99it/s]
54%|█████████████████████▋ | 543/1000 [00:53<00:41, 11.12it/s]
55%|█████████████████████▊ | 545/1000 [00:53<00:41, 10.86it/s]
55%|█████████████████████▉ | 547/1000 [00:53<00:39, 11.46it/s]
55%|█████████████████████▉ | 549/1000 [00:54<00:42, 10.63it/s]
55%|██████████████████████ | 551/1000 [00:54<00:56, 7.95it/s]
55%|██████████████████████ | 553/1000 [00:54<00:55, 8.05it/s]
55%|██████████████████████▏ | 554/1000 [00:54<00:54, 8.12it/s]
56%|██████████████████████▏ | 555/1000 [00:55<00:58, 7.61it/s]
56%|██████████████████████▏ | 556/1000 [00:55<00:55, 7.93it/s]
56%|██████████████████████▎ | 558/1000 [00:55<00:51, 8.56it/s]
56%|██████████████████████▎ | 559/1000 [00:55<00:54, 8.16it/s]
56%|██████████████████████▍ | 560/1000 [00:55<00:53, 8.24it/s]
56%|██████████████████████▍ | 561/1000 [00:55<00:53, 8.14it/s]
56%|██████████████████████▍ | 562/1000 [00:55<00:52, 8.27it/s]
56%|██████████████████████▌ | 563/1000 [00:55<00:50, 8.69it/s]
56%|██████████████████████▌ | 564/1000 [00:56<00:56, 7.74it/s]
56%|██████████████████████▌ | 565/1000 [00:56<00:59, 7.25it/s]
57%|██████████████████████▋ | 566/1000 [00:56<00:58, 7.36it/s]
57%|██████████████████████▋ | 567/1000 [00:56<00:57, 7.54it/s]
57%|██████████████████████▊ | 569/1000 [00:56<00:48, 8.94it/s]
57%|██████████████████████▊ | 571/1000 [00:56<00:43, 9.98it/s]
57%|██████████████████████▉ | 573/1000 [00:57<00:42, 9.94it/s]
57%|██████████████████████▉ | 574/1000 [00:57<00:44, 9.65it/s]
58%|███████████████████████ | 576/1000 [00:57<00:40, 10.50it/s]
58%|███████████████████████ | 578/1000 [00:57<00:38, 11.05it/s]
58%|███████████████████████▏ | 580/1000 [00:57<00:36, 11.54it/s]
58%|███████████████████████▎ | 582/1000 [00:57<00:34, 12.21it/s]
58%|███████████████████████▎ | 584/1000 [00:57<00:34, 12.14it/s]
59%|███████████████████████▍ | 586/1000 [00:58<00:32, 12.63it/s]
59%|███████████████████████▌ | 588/1000 [00:58<00:32, 12.77it/s]
59%|███████████████████████▌ | 590/1000 [00:58<00:31, 12.91it/s]
59%|███████████████████████▋ | 592/1000 [00:58<00:30, 13.20it/s]
59%|███████████████████████▊ | 594/1000 [00:58<00:31, 12.69it/s]
60%|███████████████████████▊ | 596/1000 [00:58<00:34, 11.67it/s]
60%|███████████████████████▉ | 598/1000 [00:59<00:33, 11.86it/s]
60%|████████████████████████ | 600/1000 [00:59<00:32, 12.42it/s]
60%|████████████████████████ | 602/1000 [00:59<00:31, 12.59it/s]
60%|████████████████████████▏ | 604/1000 [00:59<00:35, 11.14it/s]
61%|████████████████████████▏ | 606/1000 [00:59<00:35, 10.94it/s]
61%|████████████████████████▎ | 608/1000 [00:59<00:34, 11.47it/s]
61%|████████████████████████▍ | 610/1000 [01:00<00:37, 10.37it/s]
61%|████████████████████████▍ | 612/1000 [01:00<00:35, 11.05it/s]
61%|████████████████████████▌ | 614/1000 [01:00<00:34, 11.19it/s]
62%|████████████████████████▋ | 616/1000 [01:00<00:37, 10.13it/s]
62%|████████████████████████▋ | 618/1000 [01:01<00:40, 9.42it/s]
62%|████████████████████████▊ | 620/1000 [01:01<00:37, 10.09it/s]
62%|████████████████████████▉ | 622/1000 [01:01<00:34, 10.83it/s]
62%|████████████████████████▉ | 624/1000 [01:01<00:37, 9.94it/s]
63%|█████████████████████████ | 626/1000 [01:01<00:35, 10.64it/s]
63%|█████████████████████████ | 628/1000 [01:01<00:34, 10.85it/s]
63%|█████████████████████████▏ | 630/1000 [01:02<00:32, 11.38it/s]
63%|█████████████████████████▎ | 632/1000 [01:02<00:33, 11.04it/s]
63%|█████████████████████████▎ | 634/1000 [01:02<00:31, 11.63it/s]
64%|█████████████████████████▍ | 636/1000 [01:02<00:30, 11.87it/s]
64%|█████████████████████████▌ | 638/1000 [01:02<00:30, 11.86it/s]
64%|█████████████████████████▌ | 640/1000 [01:02<00:30, 11.91it/s]
64%|█████████████████████████▋ | 642/1000 [01:03<00:30, 11.87it/s]
64%|█████████████████████████▊ | 644/1000 [01:03<00:29, 11.95it/s]
65%|█████████████████████████▊ | 646/1000 [01:03<00:28, 12.59it/s]
65%|█████████████████████████▉ | 648/1000 [01:03<00:28, 12.17it/s]
65%|██████████████████████████ | 650/1000 [01:03<00:27, 12.70it/s]
65%|██████████████████████████ | 652/1000 [01:03<00:29, 11.66it/s]
65%|██████████████████████████▏ | 654/1000 [01:04<00:29, 11.83it/s]
66%|██████████████████████████▏ | 656/1000 [01:04<00:28, 12.24it/s]
66%|██████████████████████████▎ | 658/1000 [01:04<00:26, 12.73it/s]
66%|██████████████████████████▍ | 660/1000 [01:04<00:26, 13.06it/s]
66%|██████████████████████████▍ | 662/1000 [01:04<00:25, 13.49it/s]
66%|██████████████████████████▌ | 664/1000 [01:04<00:24, 13.46it/s]
67%|██████████████████████████▋ | 666/1000 [01:04<00:25, 13.31it/s]
67%|██████████████████████████▋ | 668/1000 [01:05<00:25, 12.92it/s]
67%|██████████████████████████▊ | 670/1000 [01:05<00:25, 12.73it/s]
67%|██████████████████████████▉ | 672/1000 [01:05<00:25, 13.07it/s]
67%|██████████████████████████▉ | 674/1000 [01:05<00:25, 13.00it/s]
68%|███████████████████████████ | 676/1000 [01:05<00:25, 12.87it/s]
68%|███████████████████████████ | 678/1000 [01:05<00:28, 11.27it/s]
68%|███████████████████████████▏ | 680/1000 [01:06<00:26, 12.04it/s]
68%|███████████████████████████▎ | 682/1000 [01:06<00:25, 12.26it/s]
68%|███████████████████████████▎ | 684/1000 [01:06<00:25, 12.49it/s]
69%|███████████████████████████▍ | 686/1000 [01:06<00:24, 12.78it/s]
69%|███████████████████████████▌ | 688/1000 [01:06<00:24, 12.84it/s]
69%|███████████████████████████▌ | 690/1000 [01:06<00:26, 11.59it/s]
69%|███████████████████████████▋ | 692/1000 [01:07<00:27, 11.05it/s]
69%|███████████████████████████▊ | 694/1000 [01:07<00:26, 11.70it/s]
70%|███████████████████████████▊ | 696/1000 [01:07<00:25, 12.11it/s]
70%|███████████████████████████▉ | 698/1000 [01:07<00:24, 12.52it/s]
70%|████████████████████████████ | 700/1000 [01:07<00:24, 12.44it/s]
70%|████████████████████████████ | 702/1000 [01:07<00:23, 12.91it/s]
70%|████████████████████████████▏ | 704/1000 [01:08<00:23, 12.84it/s]
71%|████████████████████████████▏ | 706/1000 [01:08<00:22, 13.04it/s]
71%|████████████████████████████▎ | 708/1000 [01:08<00:23, 12.68it/s]
71%|████████████████████████████▍ | 710/1000 [01:08<00:22, 13.02it/s]
71%|████████████████████████████▍ | 712/1000 [01:08<00:23, 12.40it/s]
71%|████████████████████████████▌ | 714/1000 [01:08<00:22, 12.83it/s]
72%|████████████████████████████▋ | 716/1000 [01:08<00:21, 13.13it/s]
72%|████████████████████████████▋ | 718/1000 [01:09<00:21, 13.40it/s]
72%|████████████████████████████▊ | 720/1000 [01:09<00:20, 13.49it/s]
72%|████████████████████████████▉ | 722/1000 [01:09<00:21, 12.97it/s]
72%|████████████████████████████▉ | 724/1000 [01:09<00:21, 12.87it/s]
73%|█████████████████████████████ | 726/1000 [01:09<00:21, 12.46it/s]
73%|█████████████████████████████ | 728/1000 [01:09<00:21, 12.67it/s]
73%|█████████████████████████████▏ | 730/1000 [01:10<00:20, 12.88it/s]
73%|█████████████████████████████▎ | 732/1000 [01:10<00:20, 12.78it/s]
73%|█████████████████████████████▎ | 734/1000 [01:10<00:20, 12.73it/s]
74%|█████████████████████████████▍ | 736/1000 [01:10<00:20, 12.66it/s]
74%|█████████████████████████████▌ | 738/1000 [01:10<00:20, 12.62it/s]
74%|█████████████████████████████▌ | 740/1000 [01:10<00:20, 12.69it/s]
74%|█████████████████████████████▋ | 742/1000 [01:10<00:20, 12.81it/s]
74%|█████████████████████████████▊ | 744/1000 [01:11<00:22, 11.29it/s]
75%|█████████████████████████████▊ | 746/1000 [01:11<00:21, 11.80it/s]
75%|█████████████████████████████▉ | 748/1000 [01:11<00:20, 12.38it/s]
75%|██████████████████████████████ | 750/1000 [01:11<00:21, 11.52it/s]
75%|██████████████████████████████ | 752/1000 [01:11<00:20, 11.85it/s]
75%|██████████████████████████████▏ | 754/1000 [01:12<00:20, 12.29it/s]
76%|██████████████████████████████▏ | 756/1000 [01:12<00:20, 12.03it/s]
76%|██████████████████████████████▎ | 758/1000 [01:12<00:20, 11.93it/s]
76%|██████████████████████████████▍ | 760/1000 [01:12<00:20, 12.00it/s]
76%|██████████████████████████████▍ | 762/1000 [01:12<00:18, 12.58it/s]
76%|██████████████████████████████▌ | 764/1000 [01:12<00:18, 12.68it/s]
77%|██████████████████████████████▋ | 766/1000 [01:12<00:18, 12.53it/s]
77%|██████████████████████████████▋ | 768/1000 [01:13<00:19, 11.77it/s]
77%|██████████████████████████████▊ | 770/1000 [01:13<00:20, 11.48it/s]
77%|██████████████████████████████▉ | 772/1000 [01:13<00:18, 12.19it/s]
77%|██████████████████████████████▉ | 774/1000 [01:13<00:17, 12.79it/s]
78%|███████████████████████████████ | 776/1000 [01:13<00:17, 12.82it/s]
78%|███████████████████████████████ | 778/1000 [01:13<00:17, 12.77it/s]
78%|███████████████████████████████▏ | 780/1000 [01:14<00:17, 12.71it/s]
78%|███████████████████████████████▎ | 782/1000 [01:14<00:21, 10.16it/s]
78%|███████████████████████████████▎ | 784/1000 [01:14<00:23, 9.13it/s]
78%|███████████████████████████████▍ | 785/1000 [01:14<00:25, 8.57it/s]
79%|███████████████████████████████▍ | 786/1000 [01:14<00:26, 7.95it/s]
79%|███████████████████████████████▍ | 787/1000 [01:15<00:29, 7.15it/s]
79%|███████████████████████████████▌ | 788/1000 [01:15<00:29, 7.11it/s]
79%|███████████████████████████████▌ | 789/1000 [01:15<00:29, 7.15it/s]
79%|███████████████████████████████▌ | 790/1000 [01:15<00:28, 7.43it/s]
79%|███████████████████████████████▋ | 792/1000 [01:15<00:24, 8.62it/s]
79%|███████████████████████████████▊ | 794/1000 [01:16<00:24, 8.24it/s]
80%|███████████████████████████████▊ | 795/1000 [01:16<00:25, 7.91it/s]
80%|███████████████████████████████▊ | 796/1000 [01:16<00:28, 7.25it/s]
80%|███████████████████████████████▉ | 797/1000 [01:16<00:29, 6.88it/s]
80%|███████████████████████████████▉ | 798/1000 [01:16<00:31, 6.50it/s]
80%|███████████████████████████████▉ | 799/1000 [01:16<00:31, 6.42it/s]
80%|████████████████████████████████ | 800/1000 [01:17<00:31, 6.30it/s]
80%|████████████████████████████████ | 801/1000 [01:17<00:30, 6.60it/s]
80%|████████████████████████████████ | 802/1000 [01:17<00:28, 6.87it/s]
80%|████████████████████████████████ | 803/1000 [01:17<00:28, 6.96it/s]
80%|████████████████████████████████▏ | 804/1000 [01:17<00:30, 6.51it/s]
80%|████████████████████████████████▏ | 805/1000 [01:17<00:29, 6.56it/s]
81%|████████████████████████████████▏ | 806/1000 [01:17<00:29, 6.49it/s]
81%|████████████████████████████████▎ | 807/1000 [01:18<00:29, 6.47it/s]
81%|████████████████████████████████▎ | 808/1000 [01:18<00:30, 6.30it/s]
81%|████████████████████████████████▎ | 809/1000 [01:18<00:41, 4.63it/s]
81%|████████████████████████████████▍ | 810/1000 [01:18<00:36, 5.25it/s]
81%|████████████████████████████████▍ | 812/1000 [01:18<00:27, 6.82it/s]
81%|████████████████████████████████▌ | 813/1000 [01:18<00:25, 7.34it/s]
82%|████████████████████████████████▌ | 815/1000 [01:19<00:20, 8.95it/s]
82%|████████████████████████████████▋ | 817/1000 [01:19<00:18, 9.85it/s]
82%|████████████████████████████████▊ | 819/1000 [01:19<00:17, 10.19it/s]
82%|████████████████████████████████▊ | 821/1000 [01:19<00:17, 10.05it/s]
82%|████████████████████████████████▉ | 823/1000 [01:19<00:17, 9.93it/s]
82%|█████████████████████████████████ | 825/1000 [01:20<00:16, 10.49it/s]
83%|█████████████████████████████████ | 827/1000 [01:20<00:15, 11.07it/s]
83%|█████████████████████████████████▏ | 829/1000 [01:20<00:14, 11.44it/s]
83%|█████████████████████████████████▏ | 831/1000 [01:20<00:14, 11.28it/s]
83%|█████████████████████████████████▎ | 833/1000 [01:20<00:14, 11.66it/s]
84%|█████████████████████████████████▍ | 835/1000 [01:20<00:13, 11.90it/s]
84%|█████████████████████████████████▍ | 837/1000 [01:21<00:13, 12.16it/s]
84%|█████████████████████████████████▌ | 839/1000 [01:21<00:13, 12.18it/s]
84%|█████████████████████████████████▋ | 841/1000 [01:21<00:15, 10.12it/s]
84%|█████████████████████████████████▋ | 843/1000 [01:21<00:16, 9.37it/s]
84%|█████████████████████████████████▊ | 845/1000 [01:21<00:15, 10.03it/s]
85%|█████████████████████████████████▉ | 847/1000 [01:22<00:14, 10.27it/s]
85%|█████████████████████████████████▉ | 849/1000 [01:22<00:13, 10.83it/s]
85%|██████████████████████████████████ | 851/1000 [01:22<00:13, 11.19it/s]
85%|██████████████████████████████████ | 853/1000 [01:22<00:12, 11.32it/s]
86%|██████████████████████████████████▏ | 855/1000 [01:22<00:12, 11.18it/s]
86%|██████████████████████████████████▎ | 857/1000 [01:22<00:12, 11.20it/s]
86%|██████████████████████████████████▎ | 859/1000 [01:23<00:12, 11.61it/s]
86%|██████████████████████████████████▍ | 861/1000 [01:23<00:11, 11.93it/s]
86%|██████████████████████████████████▌ | 863/1000 [01:23<00:11, 12.24it/s]
86%|██████████████████████████████████▌ | 865/1000 [01:23<00:11, 11.81it/s]
87%|██████████████████████████████████▋ | 867/1000 [01:23<00:11, 11.99it/s]
87%|██████████████████████████████████▊ | 869/1000 [01:24<00:12, 10.71it/s]
87%|██████████████████████████████████▊ | 871/1000 [01:24<00:11, 10.84it/s]
87%|██████████████████████████████████▉ | 873/1000 [01:24<00:11, 11.34it/s]
88%|███████████████████████████████████ | 875/1000 [01:24<00:11, 11.25it/s]
88%|███████████████████████████████████ | 877/1000 [01:24<00:11, 10.92it/s]
88%|███████████████████████████████████▏ | 879/1000 [01:24<00:11, 10.83it/s]
88%|███████████████████████████████████▏ | 881/1000 [01:25<00:11, 10.73it/s]
88%|███████████████████████████████████▎ | 883/1000 [01:25<00:10, 10.68it/s]
88%|███████████████████████████████████▍ | 885/1000 [01:25<00:10, 10.67it/s]
89%|███████████████████████████████████▍ | 887/1000 [01:25<00:10, 10.52it/s]
89%|███████████████████████████████████▌ | 889/1000 [01:25<00:11, 9.75it/s]
89%|███████████████████████████████████▌ | 890/1000 [01:26<00:11, 9.39it/s]
89%|███████████████████████████████████▋ | 891/1000 [01:26<00:12, 8.72it/s]
89%|███████████████████████████████████▋ | 892/1000 [01:26<00:12, 8.94it/s]
89%|███████████████████████████████████▊ | 894/1000 [01:26<00:10, 10.00it/s]
90%|███████████████████████████████████▊ | 896/1000 [01:26<00:10, 10.37it/s]
90%|███████████████████████████████████▉ | 898/1000 [01:26<00:10, 10.14it/s]
90%|████████████████████████████████████ | 900/1000 [01:27<00:09, 10.52it/s]
90%|████████████████████████████████████ | 902/1000 [01:27<00:09, 10.18it/s]
90%|████████████████████████████████████▏ | 904/1000 [01:27<00:08, 10.74it/s]
91%|████████████████████████████████████▏ | 906/1000 [01:27<00:08, 10.93it/s]
91%|████████████████████████████████████▎ | 908/1000 [01:27<00:08, 10.62it/s]
91%|████████████████████████████████████▍ | 910/1000 [01:27<00:08, 10.56it/s]
91%|████████████████████████████████████▍ | 912/1000 [01:28<00:07, 11.25it/s]
91%|████████████████████████████████████▌ | 914/1000 [01:28<00:07, 11.56it/s]
92%|████████████████████████████████████▋ | 916/1000 [01:28<00:07, 11.66it/s]
92%|████████████████████████████████████▋ | 918/1000 [01:28<00:06, 11.89it/s]
92%|████████████████████████████████████▊ | 920/1000 [01:28<00:06, 11.92it/s]
92%|████████████████████████████████████▉ | 922/1000 [01:28<00:06, 12.13it/s]
92%|████████████████████████████████████▉ | 924/1000 [01:29<00:06, 12.28it/s]
93%|█████████████████████████████████████ | 926/1000 [01:29<00:05, 12.59it/s]
93%|█████████████████████████████████████ | 928/1000 [01:29<00:05, 12.73it/s]
93%|█████████████████████████████████████▏ | 930/1000 [01:29<00:05, 12.74it/s]
93%|█████████████████████████████████████▎ | 932/1000 [01:29<00:05, 12.22it/s]
93%|█████████████████████████████████████▎ | 934/1000 [01:29<00:05, 11.96it/s]
94%|█████████████████████████████████████▍ | 936/1000 [01:30<00:05, 12.09it/s]
94%|█████████████████████████████████████▌ | 938/1000 [01:30<00:05, 10.96it/s]
94%|█████████████████████████████████████▌ | 940/1000 [01:30<00:05, 10.09it/s]
94%|█████████████████████████████████████▋ | 942/1000 [01:30<00:06, 9.44it/s]
94%|█████████████████████████████████████▊ | 944/1000 [01:30<00:05, 9.89it/s]
95%|█████████████████████████████████████▊ | 946/1000 [01:31<00:05, 10.63it/s]
95%|█████████████████████████████████████▉ | 948/1000 [01:31<00:04, 10.59it/s]
95%|██████████████████████████████████████ | 950/1000 [01:31<00:06, 8.20it/s]
95%|██████████████████████████████████████ | 952/1000 [01:31<00:05, 9.17it/s]
95%|██████████████████████████████████████▏ | 954/1000 [01:31<00:04, 9.90it/s]
96%|██████████████████████████████████████▏ | 956/1000 [01:32<00:04, 10.56it/s]
96%|██████████████████████████████████████▎ | 958/1000 [01:32<00:03, 10.65it/s]
96%|██████████████████████████████████████▍ | 960/1000 [01:32<00:04, 9.81it/s]
96%|██████████████████████████████████████▍ | 962/1000 [01:32<00:03, 10.12it/s]
96%|██████████████████████████████████████▌ | 964/1000 [01:32<00:03, 10.55it/s]
97%|██████████████████████████████████████▋ | 966/1000 [01:33<00:03, 10.86it/s]
97%|██████████████████████████████████████▋ | 968/1000 [01:33<00:02, 11.34it/s]
97%|██████████████████████████████████████▊ | 970/1000 [01:33<00:02, 11.50it/s]
97%|██████████████████████████████████████▉ | 972/1000 [01:33<00:02, 11.19it/s]
97%|██████████████████████████████████████▉ | 974/1000 [01:33<00:02, 10.37it/s]
98%|███████████████████████████████████████ | 976/1000 [01:34<00:02, 10.19it/s]
98%|███████████████████████████████████████ | 978/1000 [01:34<00:02, 10.27it/s]
98%|███████████████████████████████████████▏| 980/1000 [01:34<00:01, 10.64it/s]
98%|███████████████████████████████████████▎| 982/1000 [01:34<00:01, 10.70it/s]
98%|███████████████████████████████████████▎| 984/1000 [01:34<00:01, 11.02it/s]
99%|███████████████████████████████████████▍| 986/1000 [01:34<00:01, 10.72it/s]
99%|███████████████████████████████████████▌| 988/1000 [01:35<00:01, 10.88it/s]
99%|███████████████████████████████████████▌| 990/1000 [01:35<00:00, 11.03it/s]
99%|███████████████████████████████████████▋| 992/1000 [01:35<00:00, 11.31it/s]
99%|███████████████████████████████████████▊| 994/1000 [01:35<00:00, 11.54it/s]
100%|███████████████████████████████████████▊| 996/1000 [01:35<00:00, 11.66it/s]
100%|███████████████████████████████████████▉| 998/1000 [01:35<00:00, 11.49it/s]
100%|███████████████████████████████████████| 1000/1000 [01:36<00:00, 10.40it/s]
[54]:
import matplotlib.pyplot as plt
plt.style.use('MyStyle')
plt.clf()
plt.plot(thetas[:,0], outpar[:,0], 'o', alpha=0.4)
plt.xlabel('True $\mathcal{M}_{3D}$')
plt.ylabel('Recovered $\mathcal{M}_{3D}$')
xm = np.linspace(0., 2., 100)
plt.plot(xm, xm, ':', c='k')
[54]:
[<matplotlib.lines.Line2D at 0x1674ef690>]
[55]:
plt.clf()
plt.plot(thetas[:,1], outpar[:,1], 'o', alpha=0.4)
plt.xlabel('True $\\alpha$')
plt.ylabel('Recovered $\\alpha$')
xa = np.linspace(2., 8., 100)
plt.plot(xa, xa, ':', c='k')
[55]:
[<matplotlib.lines.Line2D at 0x167676910>]
[56]:
plt.clf()
plt.plot(thetas[:,2], outpar[:,2], 'o', alpha=0.4)
plt.xlabel('True $k_{inj}$')
plt.ylabel('Recovered $k_{inj}$')
xk = np.linspace(-4., -2., 100)
plt.plot(xk, xk, ':', c='k')
[56]:
[<matplotlib.lines.Line2D at 0x1677cc250>]
[60]:
from scipy.stats import gaussian_kde
values = np.vstack([thetas[:,0], outpar[:,0]])
kernel = gaussian_kde(values)
xmin = np.min(values[0])
xmax = np.max(values[0])
ymin = np.min(values[1])
ymax = np.max(values[1])
X, Y = np.mgrid[xmin:xmax:100j, ymin:ymax:100j]
positions = np.vstack([X.ravel(), Y.ravel()])
Z = np.reshape(kernel(positions).T, X.shape)
plt.clf()
plt.scatter(thetas[:,0], outpar[:,0], c='k', alpha=0.2, s=2)
plt.imshow(np.rot90(Z), cmap=plt.cm.gist_earth_r,
extent=[xmin, xmax, ymin, ymax], aspect='auto')
plt.xlabel('True $\mathcal{M}_{3D}$')
plt.ylabel('Recovered $\mathcal{M}_{3D}$')
plt.plot(xm, xm, ':', c='k')
[60]:
[<matplotlib.lines.Line2D at 0x16a294290>]
[61]:
values = np.vstack([thetas[:,1], outpar[:,1]])
kernel = gaussian_kde(values)
xmin = np.min(values[0])
xmax = np.max(values[0])
ymin = np.min(values[1])
ymax = np.max(values[1])
X, Y = np.mgrid[xmin:xmax:100j, ymin:ymax:100j]
positions = np.vstack([X.ravel(), Y.ravel()])
Z = np.reshape(kernel(positions).T, X.shape)
plt.clf()
plt.scatter(thetas[:,1], outpar[:,1], c='k', alpha=0.2, s=2)
plt.imshow(np.rot90(Z), cmap=plt.cm.gist_earth_r,
extent=[xmin, xmax, ymin, ymax], aspect='auto')
plt.xlabel('True $\\alpha$')
plt.ylabel('Recovered $\\alpha$')
plt.plot(xa, xa, ':', c='k')
[61]:
[<matplotlib.lines.Line2D at 0x16a41a710>]
[63]:
values = np.vstack([thetas[:,2], outpar[:,2]])
kernel = gaussian_kde(values)
xmin = np.min(values[0])
xmax = np.max(values[0])
ymin = np.min(values[1])
ymax = np.max(values[1])
X, Y = np.mgrid[xmin:xmax:100j, ymin:ymax:100j]
positions = np.vstack([X.ravel(), Y.ravel()])
Z = np.reshape(kernel(positions).T, X.shape)
plt.clf()
plt.scatter(thetas[:,2], outpar[:,2], c='k', alpha=0.2, s=2)
plt.imshow(np.rot90(Z), cmap=plt.cm.gist_earth_r,
extent=[xmin, xmax, ymin, ymax], aspect='auto')
plt.xlabel('True $k_{inj}$')
plt.ylabel('Recovered $k_{inj}$')
plt.plot(xk, xk, ':', c='k')
[63]:
[<matplotlib.lines.Line2D at 0x16a56a950>]
[100]:
med = np.median(mod.chains, axis=0)
pct = np.percentile(mod.chains, [16., 84.], axis=0)
[87]:
pct
[87]:
array([[ 0.49367922, 3.43265095, -3.83125573],
[ 1.09313597, 7.62941473, -3.03908404]])
[70]:
train_thetas = np.loadtxt('PS/sbi_model_theta.dat')
train_x = np.loadtxt('PS/sbi_model_sims.dat')
[101]:
similar = np.where((train_thetas[:,0]>=pct[0,0]) & (train_thetas[:,0]<=pct[1,0]) & (train_thetas[:,1]>=pct[0,1]) & (train_thetas[:,1]<=pct[1,1]) & (train_thetas[:,2]>=pct[0,2]) & (train_thetas[:,2]<=pct[1,2]))
[102]:
x_similar = train_x[similar]
[103]:
similar
[103]:
(array([ 56, 74, 267, 296, 328, 329, 405, 406, 475, 497, 569,
579, 600, 652, 852, 874, 892, 921, 925, 952, 962, 1072,
1080, 1108, 1145, 1147, 1221, 1243, 1253, 1257, 1313, 1332, 1350,
1380, 1388, 1447, 1487, 1492, 1546, 1555, 1605, 1735, 1819, 1872,
1961, 1965, 2031, 2037, 2075, 2111, 2246, 2249, 2263, 2315, 2341,
2378, 2429, 2474, 2475, 2480, 2495, 2497, 2511, 2517, 2614, 2654,
2712, 2736, 2839, 3007, 3025, 3097, 3123, 3214, 3287, 3291, 3440,
3452, 3460, 3470, 3514, 3538, 3594, 3749, 3846, 3862, 3883, 3891,
3899, 3949, 4095, 4100, 4147, 4205, 4242, 4250, 4262, 4267, 4377,
4378, 4393, 4451, 4551, 4640, 4722, 4724, 4836, 4861, 4889, 4919,
4945, 4964, 5056, 5075, 5089, 5150, 5163, 5189, 5258, 5324, 5351,
5360, 5410, 5426, 5471, 5567, 5625, 5753, 5778, 5806, 5864, 5916,
5973, 5999, 6048, 6255, 6310, 6470, 6477, 6494, 6499, 6513, 6555,
6612, 6619, 6632, 6655, 6685, 6757, 6827, 6862, 6922, 6948, 6969,
7014, 7017, 7042, 7089, 7143, 7196, 7232, 7292, 7296, 7307, 7368,
7390, 7418, 7441, 7451, 7453, 7456, 7512, 7534, 7573, 7595, 7607,
7657, 7679, 7739, 7829, 7834, 7874, 7898, 7906, 7963, 8201, 8264,
8307, 8336, 8436, 8452, 8578, 8631, 8677, 8689, 8699, 8894, 9024,
9182, 9252, 9366, 9383, 9405, 9477, 9488, 9536, 9563, 9619, 9660,
9678, 9715, 9799, 9915, 9919, 9963, 9985]),)
[121]:
plt.clf()
med_similar = np.median(np.abs(x_similar), axis=0)
low_similar = np.percentile(np.abs(x_similar), 16., axis=0)
high_similar = np.percentile(np.abs(x_similar), 84., axis=0)
data = np.append(np.abs(cps.bulkvel), cps.vdisp)
errs = np.append(cps.bulkerrors, cps.vdisperrors)
xr = ['$v_{b,1}$', '$v_{b,2}$', '$v_{b,3}$', '$v_{b,4}$', '$v_{b,5}$', '$v_{b,6}$', '$v_{b,7}$', '$v_{b,8}$', '$\\sigma_{v,1}$', '$\\sigma_{v,2}$', '$\\sigma_{v,3}$', '$\\sigma_{v,4}$', '$\\sigma_{v,5}$']
plt.errorbar(xr, med_similar, yerr=[med_similar-low_similar, high_similar-med_similar], fmt='s', alpha=0.5, label='Similar simulations')
plt.errorbar(xr, data, yerr=errs, fmt='o', label='Data')
plt.legend(fontsize=18)
plt.ylabel('Velocity [km/s]')
[121]:
Text(0, 0.5, 'Velocity [km/s]')
[138]:
np.savetxt('PS/tarp_data.dat', np.hstack([thetas, outpar]))
[ ]: