-
Notifications
You must be signed in to change notification settings - Fork 1
Commit
This commit does not belong to any branch on this repository, and may belong to a fork outside of the repository.
- Loading branch information
1 parent
02b2806
commit 9ea8dfd
Showing
6 changed files
with
2,436 additions
and
23 deletions.
There are no files selected for viewing
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
|
@@ -4,6 +4,7 @@ | |
*.npz | ||
*.fits | ||
*.FITS | ||
*.ascii | ||
docs/build/ | ||
|
||
!tests/data/header_only.fits | ||
|
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,316 @@ | ||
{ | ||
"cells": [ | ||
{ | ||
"cell_type": "markdown", | ||
"metadata": {}, | ||
"source": [ | ||
"# Example\n", | ||
"---" | ||
] | ||
}, | ||
{ | ||
"cell_type": "code", | ||
"execution_count": 31, | ||
"metadata": {}, | ||
"outputs": [], | ||
"source": [ | ||
"import torch\n", | ||
"import numpy as np\n", | ||
"\n", | ||
"from astropy import units, constants\n", | ||
"\n", | ||
"from p3droslo.plot import plot_cube_2D\n", | ||
"from p3droslo.model import TensorModel\n", | ||
"from p3droslo.lines import Line\n", | ||
"from p3droslo.haar import Haar" | ||
] | ||
}, | ||
{ | ||
"cell_type": "markdown", | ||
"metadata": {}, | ||
"source": [ | ||
"## Download the data\n", | ||
"---" | ||
] | ||
}, | ||
{ | ||
"cell_type": "code", | ||
"execution_count": 5, | ||
"metadata": {}, | ||
"outputs": [], | ||
"source": [ | ||
"# link to the model, and file name\n", | ||
"input_link = \"https://owncloud.ster.kuleuven.be/index.php/s/6mCZjZ2erTsXq5Y/download\"\n", | ||
"input_file = \"model_Phantom_3D.ascii\"\n", | ||
"# link to the molecular line data, and file name\n", | ||
"lamda_link = \"https://home.strw.leidenuniv.nl/~moldata/datafiles/co.dat\"\n", | ||
"lamda_file = \"co.dat\"" | ||
] | ||
}, | ||
{ | ||
"cell_type": "code", | ||
"execution_count": 6, | ||
"metadata": {}, | ||
"outputs": [ | ||
{ | ||
"name": "stdout", | ||
"output_type": "stream", | ||
"text": [ | ||
"--2023-09-27 12:41:54-- https://owncloud.ster.kuleuven.be/index.php/s/6mCZjZ2erTsXq5Y/download\n", | ||
"Resolving owncloud.ster.kuleuven.be (owncloud.ster.kuleuven.be)... 134.58.130.75\n", | ||
"Connecting to owncloud.ster.kuleuven.be (owncloud.ster.kuleuven.be)|134.58.130.75|:443... connected.\n", | ||
"HTTP request sent, awaiting response... 200 OK\n", | ||
"Length: 388714637 (371M) [application/octet-stream]\n", | ||
"Saving to: ‘model_Phantom_3D.ascii’\n", | ||
"\n", | ||
"model_Phantom_3D.as 100%[===================>] 370.71M 67.4MB/s in 5.5s \n", | ||
"\n", | ||
"2023-09-27 12:42:00 (67.0 MB/s) - ‘model_Phantom_3D.ascii’ saved [388714637/388714637]\n", | ||
"\n", | ||
"--2023-09-27 12:42:09-- https://home.strw.leidenuniv.nl/~moldata/datafiles/co.dat\n", | ||
"Resolving home.strw.leidenuniv.nl (home.strw.leidenuniv.nl)... 132.229.214.179\n", | ||
"Connecting to home.strw.leidenuniv.nl (home.strw.leidenuniv.nl)|132.229.214.179|:443... connected.\n", | ||
"HTTP request sent, awaiting response... 200 OK\n", | ||
"Length: 444204 (434K)\n", | ||
"Saving to: ‘co.dat’\n", | ||
"\n", | ||
"co.dat 100%[===================>] 433.79K 1.69MB/s in 0.3s \n", | ||
"\n", | ||
"2023-09-27 12:42:12 (1.69 MB/s) - ‘co.dat’ saved [444204/444204]\n", | ||
"\n" | ||
] | ||
} | ||
], | ||
"source": [ | ||
"# Download the model and the line data\n", | ||
"!wget $input_link --output-document $input_file\n", | ||
"!wget $lamda_link --output-document $lamda_file" | ||
] | ||
}, | ||
{ | ||
"cell_type": "code", | ||
"execution_count": 10, | ||
"metadata": {}, | ||
"outputs": [], | ||
"source": [ | ||
"# Read the Phantom ascii file\n", | ||
"(x,y,z, h, rho, vx,vy,vz, u) = np.loadtxt(input_file, skiprows=14, usecols=(0,1,2,4,5,6,7,8,9), unpack=True)\n", | ||
"\n", | ||
"# Constants that can be read from ascii file\n", | ||
"velocity_cte = 2.9784608e+06\n", | ||
"density_cte = 5.9410314e-07\n", | ||
"energy_cte = 8.8712277e+12\n", | ||
"\n", | ||
"keep = np.logical_and(h>0.0, rho>0.0)\n", | ||
"\n", | ||
"# Exclude unphysical points and points with zero abundance\n", | ||
"x = x [keep]\n", | ||
"y = y [keep]\n", | ||
"z = z [keep]\n", | ||
"vx = vx [keep]\n", | ||
"vy = vy [keep]\n", | ||
"vz = vz [keep]\n", | ||
"u = u [keep]\n", | ||
"rho = rho[keep]\n", | ||
"\n", | ||
"# Extract the number of points\n", | ||
"npoints = len(x)\n", | ||
"\n", | ||
"# Convert rho (total density) to abundances\n", | ||
"nH2 = rho * density_cte * 1.0e+6 * constants.N_A.si.value / 2.02\n", | ||
"nCO = nH2 * 1.0e-4\n", | ||
"\n", | ||
"position = np.array((x, y, z )).T\n", | ||
"velocity = np.array((vx,vy,vz)).T\n", | ||
"\n", | ||
"# Convert units\n", | ||
"position *= constants.au.si.value # convert au to m\n", | ||
"velocity *= velocity_cte * 1.0e-2 # convert cm/s to m/s\n", | ||
"\n", | ||
"# Derive temperature from internal energy (assuming adiabatic heating/cooling)\n", | ||
"gamma = 1.2\n", | ||
"mu = 2.381\n", | ||
"tmp = mu * (gamma-1.0) * u * energy_cte * 1.00784 * (units.erg/units.g * constants.u/constants.k_B).to(units.K).value\n", | ||
"\n", | ||
"# Clamp temperatures below 2.725 K\n", | ||
"tmp[tmp<2.725] = 2.725\n", | ||
"\n", | ||
"# Define turbulence at 150 m/s\n", | ||
"trb = 150.0 * np.ones(npoints)" | ||
] | ||
}, | ||
{ | ||
"cell_type": "code", | ||
"execution_count": 34, | ||
"metadata": {}, | ||
"outputs": [], | ||
"source": [ | ||
"# Convert point data to cartesian grid and put it ina TensorModel\n", | ||
"haar = Haar(position, q=8)\n", | ||
"model = TensorModel(shape=tuple([2**(haar.q-1)] * 3), sizes=haar.xyz_L)\n", | ||
"model['density' ] = haar.map_data(nCO, interpolate=True)[-1]\n", | ||
"model['temperature'] = haar.map_data(tmp, interpolate=True)[-1]\n", | ||
"model['v_turbulence'] = haar.map_data(trb, interpolate=True)[-1]\n", | ||
"model['velocity_los'] = haar.map_data(velocity[:,2], interpolate=True)[-1] # velocity along z-axis\n" | ||
] | ||
}, | ||
{ | ||
"cell_type": "code", | ||
"execution_count": 35, | ||
"metadata": {}, | ||
"outputs": [ | ||
{ | ||
"name": "stdout", | ||
"output_type": "stream", | ||
"text": [ | ||
"You have selected line:\n", | ||
" co(J=1-0)\n", | ||
"Please check the properties that were inferred:\n", | ||
" Frequency 1.152712018e+11 Hz\n", | ||
" Einstein A coeff 7.203000000e-08 1/s\n", | ||
" Molar mass 28.0 g/mol\n" | ||
] | ||
} | ||
], | ||
"source": [ | ||
"line = Line(\n", | ||
" species_name = \"co\",\n", | ||
" transition = 0,\n", | ||
" datafile = lamda_file,\n", | ||
" molar_mass = 28.0\n", | ||
")" | ||
] | ||
}, | ||
{ | ||
"cell_type": "code", | ||
"execution_count": 36, | ||
"metadata": {}, | ||
"outputs": [], | ||
"source": [ | ||
"vdiff = 1500 # velocity increment size [m/s]\n", | ||
"nfreq = 31 # number of frequencies\n", | ||
"dd = vdiff / constants.c.si.value * nfreq\n", | ||
"fmin = line.frequency - line.frequency*dd\n", | ||
"fmax = line.frequency + line.frequency*dd\n", | ||
"\n", | ||
"frequencies = torch.linspace(fmin, fmax, nfreq)" | ||
] | ||
}, | ||
{ | ||
"cell_type": "code", | ||
"execution_count": 43, | ||
"metadata": {}, | ||
"outputs": [], | ||
"source": [ | ||
"img = line.LTE_image_along_last_axis(\n", | ||
" density = model['density' ],\n", | ||
" temperature = model['temperature' ],\n", | ||
" v_turbulence = model['v_turbulence'],\n", | ||
" velocity_los = model['velocity_los'],\n", | ||
" frequencies = frequencies,\n", | ||
" dx = model.dx(3-1)\n", | ||
")\n", | ||
"\n", | ||
"# Avoid negative values (should probably avoid these earlier...)\n", | ||
"img = torch.abs(img)" | ||
] | ||
}, | ||
{ | ||
"cell_type": "code", | ||
"execution_count": 46, | ||
"metadata": {}, | ||
"outputs": [ | ||
{ | ||
"data": { | ||
"application/vnd.jupyter.widget-view+json": { | ||
"model_id": "be3aff1bc931448daa764e6346cc0d0f", | ||
"version_major": 2, | ||
"version_minor": 0 | ||
}, | ||
"text/plain": [ | ||
"interactive(children=(IntSlider(value=15, description='z', max=30), Output()), _dom_classes=('widget-interact'…" | ||
] | ||
}, | ||
"metadata": {}, | ||
"output_type": "display_data" | ||
}, | ||
{ | ||
"data": { | ||
"text/plain": [ | ||
"<function p3droslo.plot.plot_cube_2D.<locals>.plot(z)>" | ||
] | ||
}, | ||
"execution_count": 46, | ||
"metadata": {}, | ||
"output_type": "execute_result" | ||
} | ||
], | ||
"source": [ | ||
"plot_cube_2D(np.log10(img+1.0e-20))" | ||
] | ||
}, | ||
{ | ||
"cell_type": "code", | ||
"execution_count": 49, | ||
"metadata": {}, | ||
"outputs": [ | ||
{ | ||
"data": { | ||
"application/vnd.jupyter.widget-view+json": { | ||
"model_id": "a34b237d0f0d41a497cb2cc726b365d3", | ||
"version_major": 2, | ||
"version_minor": 0 | ||
}, | ||
"text/plain": [ | ||
"interactive(children=(IntSlider(value=63, description='i', max=127), IntSlider(value=63, description='j', max=…" | ||
] | ||
}, | ||
"metadata": {}, | ||
"output_type": "display_data" | ||
}, | ||
{ | ||
"data": { | ||
"text/plain": [ | ||
"<function p3droslo.plot.plot_spectrum.<locals>.plot(i, j)>" | ||
] | ||
}, | ||
"execution_count": 49, | ||
"metadata": {}, | ||
"output_type": "execute_result" | ||
} | ||
], | ||
"source": [ | ||
"plot_spectrum(img)" | ||
] | ||
}, | ||
{ | ||
"cell_type": "code", | ||
"execution_count": null, | ||
"metadata": {}, | ||
"outputs": [], | ||
"source": [] | ||
} | ||
], | ||
"metadata": { | ||
"kernelspec": { | ||
"display_name": "magritte", | ||
"language": "python", | ||
"name": "python3" | ||
}, | ||
"language_info": { | ||
"codemirror_mode": { | ||
"name": "ipython", | ||
"version": 3 | ||
}, | ||
"file_extension": ".py", | ||
"mimetype": "text/x-python", | ||
"name": "python", | ||
"nbconvert_exporter": "python", | ||
"pygments_lexer": "ipython3", | ||
"version": "3.9.13" | ||
}, | ||
"orig_nbformat": 4 | ||
}, | ||
"nbformat": 4, | ||
"nbformat_minor": 2 | ||
} |
Oops, something went wrong.