Brewster_v2_tutorial.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "20481507",
   "metadata": {},
   "source": [
    "# Updates Summary\n",
    "\n",
    "## 1. Introduced Class and Dictionary into Brewster\n",
    "- introducing classes for better organization and management of code.\n",
    "- Utilized dictionaries for efficient data storage and retrieval, making it easier to index parameters and add new parameters.\n",
    "\n",
    "## 2. Merged Multinest and MCMC\n",
    "- Combined the Multinest and MCMC into a unified framework.\n",
    "\n",
    "## 3. Simplified Cloud Initialization and Retrieval Setup\n",
    "- Streamlined the process for initializing cloud parameters.\n",
    "- Reduced complexity in the setup process, automatically generate p0 in MCMC, making it more user-friendly and less error-prone.\n",
    "\n",
    "## 4. Added Non-uniform Gas Parameterization\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "15eeebd5",
   "metadata": {},
   "outputs": [],
   "source": [
    "import os \n",
    "os.chdir(\"/Users/fw23aao/brewsteretrieval/brewster\")\n",
    "import utils\n",
    "import test_module\n",
    "from collections import namedtuple\n",
    "import numpy as np\n",
    "import retrieval_run\n",
    "import settings"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "063f8cbe",
   "metadata": {},
   "source": [
    "# 1. mcmc template"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "aa21a424",
   "metadata": {},
   "outputs": [],
   "source": [
    "fwhm=700\n",
    "wavelength_range=[1,2.8]\n",
    "ndata=1\n",
    "\n",
    "chemeq=0\n",
    "gaslist = ['h2o','co','tio','vo','crh','feh','na','k']\n",
    "# gaslist = ['h2o','co','tio','vo','crh','feh','k','na']\n",
    "gastype_list=['U','U','U','U','U','N','U','U']\n",
    "\n",
    "ptype=77  \n",
    "\n",
    "do_clouds=1\n",
    "npatches=2\n",
    "\n",
    "# cloudname = ['power law cloud slab','powerlaw cloud deck']  \n",
    "# cloudpacth_index=[[1],[1,2]] \n",
    "# particle_dis=['hansan','hansan']\n",
    "\n",
    "cloudname = ['power law cloud slab']  \n",
    "# cloudname = ['power law cloud slab']  #'Mie scattering cloud deck--Al2O3'\n",
    "cloudpacth_index=[[1]] \n",
    "particle_dis=['hansan']\n",
    "\n",
    "do_fudge=1\n",
    "samplemode='mcmc'\n",
    "\n",
    "instrument_instance = utils.Instrument(fwhm,wavelength_range,ndata)\n",
    "re_params = utils.Retrieval_params(samplemode,chemeq,gaslist,gastype_list,fwhm,do_fudge,ptype,do_clouds,npatches,cloudname,cloudpacth_index,particle_dis)\n",
    "model_config_instance = utils.ModelConfig(samplemode)\n",
    "io_config_instance = utils.IOConfig()\n",
    "\n",
    "\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "e3b893dd",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Instrument: \n",
      "------------\n",
      "- fwhm : 700\n",
      "- wavelength_range : [1, 2.8]\n",
      "- ndata : 1 \n",
      "- wavpoints : None \n",
      "\n"
     ]
    }
   ],
   "source": [
    "print(instrument_instance)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "8c101a13",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "ModelConfig: \n",
      "------------\n",
      "- use_disort : 0\n",
      "- do_fudge : 1\n",
      "- malk : 0\n",
      "- mch4 : 0\n",
      "- do_bff : 1\n",
      "- fresh : 0\n",
      "- xpath : ../Linelists/\n",
      "- xlist : gaslistRox.dat\n",
      "- dist : None\n",
      "- pfile : LSR1835_eqpt.dat\n",
      "\n",
      "(MCMC): \n",
      "------------\n",
      "- ndim : None\n",
      "- nwalkers : None\n",
      "- nburn : 10000\n",
      "- niter : 30000\n",
      "\n"
     ]
    }
   ],
   "source": [
    "print(model_config_instance)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "89b52da5",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "IOConfig: \n",
      "------------\n",
      "- runname : retrieval\n",
      "- outdir : N/A\n",
      "- finalout : retrieval.pk1\n",
      "- chaindump : retrieval_last_nc.pic\n",
      "- picdump : retrieval_snapshot.pic\n",
      "- statfile : status_ball_retrieval.txt\n",
      "- rfile : runtimes_retrieval.dat\n",
      "- runtest : 1\n",
      "- make_arg_pickle : 2\n",
      "\n"
     ]
    }
   ],
   "source": [
    "print(io_config_instance)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "444a4a98",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "retrieval_param: \n",
      "------------\n",
      "- gas : ['h2o', 'co', 'tio', 'vo', 'crh', 'feh', 'na', 'k']\n",
      "- refinement_params : ['logg', 'r2d2', 'dlambda', 'tolerance_parameter_1']\n",
      "- pt : ['gamma', 'Tint', 'alpha', 'lndelta', 'T1', 'T2', 'T3']\n",
      "- cloud : ['fcld', 'patch 1']\n",
      "  -- patch 1: ['tau_pcs', 'logp_pcs', 'dp_pcs', 'omega_pcs', 'alpha_pcs']\n",
      "- added_params: []\n",
      "\n"
     ]
    }
   ],
   "source": [
    "print(re_params)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "3c94ba94",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "['customized', 0]"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "re_params.dictionary['refinement_params']['params']['tolerance_parameter_1']['distribution']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "dca5f62d",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "{'gas': {'h2o': {'gastype': 'U',\n",
       "   'params': {'log_abund': {'initialization': None,\n",
       "     'distribution': ['normal', -4.0, 0.5],\n",
       "     'range': None,\n",
       "     'prior': None}}},\n",
       "  'co': {'gastype': 'U',\n",
       "   'params': {'log_abund': {'initialization': None,\n",
       "     'distribution': ['normal', -4.0, 0.5],\n",
       "     'range': None,\n",
       "     'prior': None}}},\n",
       "  'tio': {'gastype': 'U',\n",
       "   'params': {'log_abund': {'initialization': None,\n",
       "     'distribution': ['normal', -4.0, 0.5],\n",
       "     'range': None,\n",
       "     'prior': None}}},\n",
       "  'vo': {'gastype': 'U',\n",
       "   'params': {'log_abund': {'initialization': None,\n",
       "     'distribution': ['normal', -4.0, 0.5],\n",
       "     'range': None,\n",
       "     'prior': None}}},\n",
       "  'crh': {'gastype': 'U',\n",
       "   'params': {'log_abund': {'initialization': None,\n",
       "     'distribution': ['normal', -4.0, 0.5],\n",
       "     'range': None,\n",
       "     'prior': None}}},\n",
       "  'feh': {'gastype': 'N',\n",
       "   'params': {'log_abund': {'initialization': None,\n",
       "     'distribution': ['normal', -4.0, 0.5],\n",
       "     'range': None,\n",
       "     'prior': None},\n",
       "    'p_ref': {'initialization': None,\n",
       "     'distribution': ['normal', -1, 0.2],\n",
       "     'range': None,\n",
       "     'prior': None},\n",
       "    'alpha': {'initialization': None,\n",
       "     'distribution': ['uniform', 0, 1],\n",
       "     'range': None,\n",
       "     'prior': None}}},\n",
       "  'na': {'gastype': 'U',\n",
       "   'params': {'log_abund': {'initialization': None,\n",
       "     'distribution': ['normal', -4.0, 0.5],\n",
       "     'range': None,\n",
       "     'prior': None}}},\n",
       "  'k': {'gastype': 'U',\n",
       "   'params': {'log_abund': {'initialization': None,\n",
       "     'distribution': ['normal', -4.0, 0.5],\n",
       "     'range': None,\n",
       "     'prior': None}}}},\n",
       " 'refinement_params': {'params': {'logg': {'initialization': None,\n",
       "    'distribution': ['normal', 4.5, 0.1],\n",
       "    'prior': None},\n",
       "   'r2d2': {'initialization': None,\n",
       "    'distribution': ['normal', 0, 1],\n",
       "    'prior': None},\n",
       "   'dlambda': {'initialization': None,\n",
       "    'distribution': ['normal', 0, 0.001],\n",
       "    'prior': None},\n",
       "   'tolerance_parameter_1': {'initialization': None,\n",
       "    'distribution': ['customized',\n",
       "     <function utils.tolerance_parameter_customized_distribution(x)>],\n",
       "    'range': None,\n",
       "    'prior': None}}},\n",
       " 'pt': {'ptype': 77,\n",
       "  'params': {'gamma': {'initialization': None,\n",
       "    'distribution': ['normal', 50, 1],\n",
       "    'range': None,\n",
       "    'prior': None},\n",
       "   'Tint': {'initialization': None,\n",
       "    'distribution': ['normal', 1200, 200],\n",
       "    'range': None,\n",
       "    'prior': None},\n",
       "   'alpha': {'initialization': None,\n",
       "    'distribution': ['uniform', 1, 2],\n",
       "    'range': None,\n",
       "    'prior': None},\n",
       "   'lndelta': {'initialization': None,\n",
       "    'distribution': ['normal', 0, 1],\n",
       "    'range': None,\n",
       "    'prior': None},\n",
       "   'T1': {'initialization': None,\n",
       "    'distribution': ['normal', 1200, 200],\n",
       "    'range': None,\n",
       "    'prior': None},\n",
       "   'T2': {'initialization': None,\n",
       "    'distribution': ['normal', 1200, 200],\n",
       "    'range': None,\n",
       "    'prior': None},\n",
       "   'T3': {'initialization': None,\n",
       "    'distribution': ['normal', 1200, 200],\n",
       "    'range': None,\n",
       "    'prior': None}}},\n",
       " 'cloud': {'fcld': {'initialization': None,\n",
       "   'distribution': ['uniform', 0, 1],\n",
       "   'range': None,\n",
       "   'prior': None},\n",
       "  'patch 1': {'power law cloud slab': {'cloudnum': 89,\n",
       "    'cloudtype': 1,\n",
       "    'params': {'tau_pcs': {'initialization': None,\n",
       "      'distribution': ['normal', 10, 1],\n",
       "      'range': None,\n",
       "      'prior': None},\n",
       "     'logp_pcs': {'initialization': None,\n",
       "      'distribution': ['normal', -0.2, 0.5],\n",
       "      'range': None,\n",
       "      'prior': None},\n",
       "     'dp_pcs': {'initialization': None,\n",
       "      'distribution': ['customized',\n",
       "       <function utils.dp_customized_distribution(x)>],\n",
       "      'range': None,\n",
       "      'prior': None},\n",
       "     'omega_pcs': {'initialization': None,\n",
       "      'distribution': ['uniform', 0, 1],\n",
       "      'range': None,\n",
       "      'prior': None},\n",
       "     'alpha_pcs': {'initialization': None,\n",
       "      'distribution': ['normal', 0, 1],\n",
       "      'range': None,\n",
       "      'prior': None}}}}}}"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "re_params.dictionary"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "5e84503e",
   "metadata": {},
   "outputs": [],
   "source": [
    "# re_params.dictionary\n",
    "all_params,all_params_values =utils.get_all_parametres(re_params.dictionary) \n",
    "# all_params"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "52064cf2",
   "metadata": {},
   "source": [
    "### add new parameters"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "e8c0c1ed",
   "metadata": {},
   "outputs": [],
   "source": [
    "# re_params.added_params(\"new1\")\n",
    "# print(re_params)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b89fc754",
   "metadata": {},
   "source": [
    "## 1.1 initialize the parameters"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "2d6b62c8",
   "metadata": {},
   "outputs": [],
   "source": [
    "io_config_instance.outdir=\"/beegfs/car/fei/lsr1835/\"\n",
    "io_config_instance.runname='new_test'\n",
    "io_config_instance.update_dictionary()\n",
    "\n",
    "\n",
    "model_config_instance.dist= 5.689 \n",
    "model_config_instance.update_dictionary()\n",
    "\n",
    "obspec = np.asfortranarray(np.loadtxt(\"LSR1835_data_realcalib_new_trimmed.dat\",dtype='d',unpack='true'))\n",
    "\n",
    "\n",
    "args_instance = utils.ArgsGen(re_params,model_config_instance,instrument_instance,obspec)\n",
    "\n",
    "settings.init(args_instance)\n",
    "args_instance=settings.runargs\n",
    "\n",
    "# settings.init()\n",
    "# settings.runargs=utils.args_gen(re_params,model_config_instance,instrument_instance,obspec)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7f3bfd03",
   "metadata": {},
   "source": [
    "## 1.2  set initialization distribution"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 46,
   "id": "1db79d5e",
   "metadata": {},
   "outputs": [],
   "source": [
    "# r2d2 = (71492e3)**2. / (model_config_instance.dist * 3.086e+16)**2.\n",
    "# re_params.dictionary['refinement_params']['params']['r2d2']['distribution']=['normal', r2d2, 0.1*r2d2]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "id": "60166374",
   "metadata": {},
   "outputs": [],
   "source": [
    "# model_config_instance.ndim=len(all_params)\n",
    "# model_config_instance.nwalkers=len(all_params)*16\n",
    "# p0=utils.MC_P0_gen(re_params.dictionary,model_config_instance)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f236fcfc",
   "metadata": {},
   "source": [
    "## 1.3 run retrieval"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "id": "5d0d7b60",
   "metadata": {},
   "outputs": [],
   "source": [
    "# retrieval_run.brewster_reterieval_run(re_params,model_config_instance,io_config_instance)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8d4b50d1",
   "metadata": {},
   "source": [
    "# 2. multinest template"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "id": "b2392de5",
   "metadata": {},
   "outputs": [],
   "source": [
    "fwhm=700\n",
    "wavelength_range=[1,2.8]\n",
    "ndata=1\n",
    "\n",
    "\n",
    "chemeq=0\n",
    "gaslist = ['h2o','co','tio','vo','crh','feh','na','k']\n",
    "gastype_list=['U','U','U','U','U','N','U','U']\n",
    "\n",
    "ptype=77  \n",
    "do_clouds=1\n",
    "npatches=2\n",
    "cloudname = ['power law cloud slab','powerlaw cloud deck']  \n",
    "cloudpacth_index=[[1],[1,2]] \n",
    "# particle_dis=['hansan','log_normal']\n",
    "\n",
    "\n",
    "do_fudge=1\n",
    "samplemode='multinest'\n",
    "\n",
    "instrument_instance = utils.Instrument(fwhm,wavelength_range,ndata)\n",
    "re_params = utils.Retrieval_params(samplemode,chemeq,gaslist,gastype_list,fwhm,do_fudge,ptype,do_clouds,npatches,cloudname,cloudpacth_index)\n",
    "model_config_instance = utils.ModelConfig(samplemode)\n",
    "io_config_instance = utils.IOConfig()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "id": "4c2a8ae6",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Instrument: \n",
      "------------\n",
      "- fwhm : 700\n",
      "- wavelength_range : [1, 2.8]\n",
      "- ndata : 1 \n",
      "- wavpoints : None \n",
      "\n"
     ]
    }
   ],
   "source": [
    "print(instrument_instance)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "id": "b711c2c4",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "ModelConfig: \n",
      "------------\n",
      "- use_disort : 0\n",
      "- do_fudge : 1\n",
      "- malk : 0\n",
      "- mch4 : 0\n",
      "- do_bff : 1\n",
      "- fresh : 0\n",
      "- xpath : ../Linelists/\n",
      "- xlist : gaslistRox.dat\n",
      "- dist : None\n",
      "- pfile : LSR1835_eqpt.dat\n",
      "\n",
      "PyMultiNest: \n",
      "----------------------\n",
      "- LogLikelihood: None \n",
      "- Prior: None \n",
      "- ndim: None \n",
      "- nparam: None\n",
      "- n_clustering_params: None\n",
      "- wrapped_params: None\n",
      "- importance_nested_sampling: True\n",
      "- multimodal: True\n",
      "- const_efficiency_mode: False\n",
      "- n_live_points: 400\n",
      "- evidence_tolerance: 0.5\n",
      "- sampling_efficiency: 0.8\n",
      "- n_iter_before_update: 100\n",
      "- null_log_evidence: -1e+90\n",
      "- max_modes: 100\n",
      "- mode_tolerance: -1e+90\n",
      "- outputfiles_basename: \n",
      "- seed: -1\n",
      "- verbose: True\n",
      "- resume: True\n",
      "- context: 0\n",
      "- write_output: True\n",
      "- log_zero: -1e+100\n",
      "- max_iter: 0\n",
      "- init_MPI: False\n",
      "- dump_callback: None\n",
      "\n"
     ]
    }
   ],
   "source": [
    "print(model_config_instance)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "id": "9bb8d0b2",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "IOConfig: \n",
      "------------\n",
      "- runname : retrieval\n",
      "- outdir : N/A\n",
      "- finalout : retrieval.pk1\n",
      "- chaindump : retrieval_last_nc.pic\n",
      "- picdump : retrieval_snapshot.pic\n",
      "- statfile : status_ball_retrieval.txt\n",
      "- rfile : runtimes_retrieval.dat\n",
      "- runtest : 1\n",
      "- make_arg_pickle : 2\n",
      "\n"
     ]
    }
   ],
   "source": [
    "print(io_config_instance)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "id": "a01db9ca",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "retrieval_param: \n",
      "------------\n",
      "- gas : ['h2o', 'co', 'tio', 'vo', 'crh', 'feh', 'na', 'k']\n",
      "- refinement_params : ['M', 'R', 'dlambda', 'tolerance_parameter_1']\n",
      "- pt : ['gamma', 'Tint', 'alpha', 'lndelta', 'T1', 'T2', 'T3']\n",
      "- cloud : ['fcld', 'patch 1', 'patch 2']\n",
      "  -- patch 1: ['tau_pcs', 'logp_pcs', 'dp_pcs', 'omega_pcs', 'alpha_pcs', 'logp_pcd', 'dp_pcd', 'omega_pcd', 'alpha_pcd']\n",
      "  -- patch 2: ['logp_pcd', 'dp_pcd', 'omega_pcd', 'alpha_pcd']\n",
      "- added_params: []\n",
      "\n"
     ]
    }
   ],
   "source": [
    "print(re_params)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "90a6a9b1",
   "metadata": {},
   "source": [
    "## 2.1 initialize the parameters"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "id": "b1e74b33",
   "metadata": {},
   "outputs": [],
   "source": [
    "io_config_instance.outdir=\"/beegfs/car/fei/lsr1835/test/\"\n",
    "io_config_instance.runname='new_test'\n",
    "io_config_instance.update_dictionary()\n",
    "\n",
    "\n",
    "model_config_instance.const_efficiency_mode=True\n",
    "model_config_instance.sampling_efficiency=0.3\n",
    "model_config_instance.multimodal = False\n",
    "model_config_instance.log_zero= -1e90\n",
    "model_config_instance.importance_nested_sampling= False\n",
    "model_config_instance.evidence_tolerance=0.1\n",
    "model_config_instance.dist= 5.689 \n",
    "model_config_instance.update_dictionary()\n",
    "\n",
    "obspec = np.asfortranarray(np.loadtxt(\"LSR1835_data_realcalib_new_trimmed.dat\",dtype='d',unpack='true'))\n",
    "settings.init()\n",
    "settings.runargs=utils.args_gen(re_params,model_config_instance,instrument_instance,obspec)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "22dc7e95",
   "metadata": {},
   "source": [
    "## 2.2 run retrieval"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "id": "c42e170a",
   "metadata": {},
   "outputs": [],
   "source": [
    "# retrieval_run.brewster_reterieval_run(re_params,model_config_instance,io_config_instance)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "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.8.18"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}