{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Example of asynchronous requests (v > 1.1)\n",
    "- The scope of this example is to show how to request several products together so that internal resource usage is maximized\n",
    "- We extract the spectrum of the Crab in groups of 'nscw' science windows for each year from 'start_year' to 'stop_year' included\n",
    "- We use a token provided by the web interface to receive dedicated emails\n",
    "- We optionally show how to fit the spectra with a broken power law using xspec"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": [
     "parameters"
    ]
   },
   "outputs": [],
   "source": [
    "#A few input parameters\n",
    "osa_version=\"OSA11.2\"\n",
    "source_name=\"Crab\"\n",
    "nscw=10\n",
    "start_year=2017\n",
    "end_year=2018\n",
    "systematic_fraction = 0.01\n",
    "token=''\n",
    "integral_data_rights = \"public\""
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Token authentication\n",
    "\n",
    "- You can provide a valid token as explained in the 'Authentication' example or skip the following cell and continue anonymously"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# import getpass\n",
    "# token = getpass.getpass('Insert the token')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# To know details of the token\n",
    "# import oda_api.token \n",
    "# oda_api.token.decode_oda_token(token)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#We hardcode a catalog for the Crab\n",
    "api_cat={\n",
    "    \"cat_frame\": \"fk5\", \n",
    "    \"cat_coord_units\": \"deg\", \n",
    "    \"cat_column_list\": [\n",
    "        [0, 7], \n",
    "        [\"1A 0535+262\", \"Crab\"], \n",
    "        [125.4826889038086, 1358.7255859375], \n",
    "        [84.72280883789062, 83.63166809082031], \n",
    "        [26.312734603881836, 22.016284942626953], \n",
    "        [-32768, -32768], \n",
    "        [2, 2], \n",
    "        [0, 0], \n",
    "        [0.0002800000074785203, 0.0002800000074785203]], \n",
    "    \"cat_column_names\": [\n",
    "        \"meta_ID\", \n",
    "        \"src_names\", \n",
    "        \"significance\", \n",
    "        \"ra\", \n",
    "        \"dec\", \n",
    "        \"NEW_SOURCE\", \n",
    "        \"ISGRI_FLAG\", \n",
    "        \"FLAG\", \n",
    "        \"ERR_RAD\"\n",
    "    ], \n",
    "    \"cat_column_descr\": \n",
    "        [\n",
    "            [\"meta_ID\", \"<i8\"], \n",
    "            [\"src_names\", \"<U11\"], \n",
    "            [\"significance\", \"<f8\"], \n",
    "            [\"ra\", \"<f8\"], \n",
    "            [\"dec\", \"<f8\"], \n",
    "            [\"NEW_SOURCE\", \"<i8\"], \n",
    "            [\"ISGRI_FLAG\", \"<i8\"], \n",
    "            [\"FLAG\", \"<i8\"], \n",
    "            [\"ERR_RAD\", \"<f8\"]\n",
    "        ], \n",
    "    \"cat_lat_name\": \"dec\", \n",
    "    \"cat_lon_name\": \"ra\"\n",
    "}\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Let's get some logging\n",
    "\n",
    "This is to help visualizing the progress.\n",
    "\n",
    "* WARNING is the default level\n",
    "* INFO writes some more information\n",
    "* DEBUG is maily for developers and issue tracking"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import logging\n",
    "#default\n",
    "#logging.getLogger().setLevel(logging.WARNING)\n",
    "#slightly more verbose\n",
    "logging.getLogger().setLevel(logging.INFO)\n",
    "#all messages\n",
    "#logging.getLogger().setLevel(logging.DEBUG)\n",
    "\n",
    "logging.getLogger('oda_api').addHandler(logging.StreamHandler()) "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Connection to the dispatcher\n",
    "\n",
    "We build the dispatcher object, which we will use for issuing our requests to the platform. The token will be automatically discovered, if available."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "\n",
    "import numpy as np\n",
    "import json\n",
    "import oda_api.api\n",
    "\n",
    "import oda_api\n",
    "\n",
    "disp = oda_api.api.DispatcherAPI()\n",
    "disp.get_instrument_description(\"isgri\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "- Here, we collect and spectra for each year in a random sample of nscw=10 science windows\n",
    "- We use the hard-coded catalog.\n",
    "- note that we make a loop and submit the jobs without waiting for their completion\n",
    "- at each loop, we test if they completed and we count how many have finished\n",
    "- we continue to poll the dispatcher for unfinished jobs and we terminate the loop when all are done\n",
    "- In this way, we let the platform optimize our requests\n",
    "- There will be a convenience function in future versions of oda_api for this purpose"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "spectrum_results=[]\n",
    "\n",
    "disp_by_ys = {}\n",
    "data_by_ys = {}\n",
    "\n",
    "par_dict = {'RA': 83.5,\n",
    "            'DEC': 22.0,\n",
    "            'radius': 10,\n",
    "            'instrument':'isgri',\n",
    "            'product': 'isgri_spectrum',\n",
    "            'osa_version' : osa_version,\n",
    "            'product_type': 'Real',\n",
    "            'max_pointings': nscw,\n",
    "            'selected_catalog' : json.dumps(api_cat),\n",
    "            'integral_data_rights' : integral_data_rights}\n",
    "\n",
    "# Should you need to access private data, just add this option\n",
    "#,\"integral_data_rights\": \"all-private\"}\n",
    "\n",
    "if token != '':\n",
    "    par_dict.update({'token': token})\n",
    "elif disp.token is not None:\n",
    "    par_dict.update({'token': disp.token})\n",
    "\n",
    "while True:\n",
    "    spectrum_results=[]\n",
    "\n",
    "    for year in range(start_year, end_year+1): \n",
    "        T1_utc='%4d-01-01T00:00:00.0'%year\n",
    "        T2_utc='%4d-12-31T23:59:59.0'%year \n",
    "        \n",
    "        print(T1_utc,'-',T2_utc)\n",
    "\n",
    "        par_dict.update({'T1': T1_utc,\n",
    "                        'T2': T2_utc})\n",
    "        \n",
    "        if year >= 2016:\n",
    "            osa_version='OSA11.2'\n",
    "        else:\n",
    "            osa_version='OSA10.2'\n",
    "\n",
    "        #Just renaiming for a general dictionary key\n",
    "        ys = year\n",
    "\n",
    "        # We start one dipatcher for each job,\n",
    "        # they will run in parallel until products are ready\n",
    "        if ys not in disp_by_ys:\n",
    "            disp_by_ys[ys] = oda_api.api.DispatcherAPI(url=disp.url, wait=False) #Note the flag wait=False\n",
    "\n",
    "        _disp = disp_by_ys[ys]\n",
    "\n",
    "        data = data_by_ys.get(ys, None)\n",
    "\n",
    "        if data is None and not _disp.is_failed:\n",
    "            \n",
    "            #We submit or we poll \n",
    "            if not _disp.is_submitted:\n",
    "                data = _disp.get_product(**par_dict)\n",
    "            else:\n",
    "                _disp.poll()\n",
    "\n",
    "            print(\"Is complete \", _disp.is_complete)\n",
    "            # We retrieve data\n",
    "            if not _disp.is_complete:\n",
    "                continue\n",
    "            else:\n",
    "                data = _disp.get_product(**par_dict)\n",
    "                data_by_ys[ys] = data                \n",
    "\n",
    "        spectrum_results.append(data)\n",
    "        \n",
    "    n_complete = len([ year for year, _disp in disp_by_ys.items() if _disp.is_complete ])\n",
    "    print(f\"complete {n_complete} / {len(disp_by_ys)}\")\n",
    "    if n_complete == len(disp_by_ys):\n",
    "        print(\"done!\")\n",
    "        break\n",
    "    print(\"not done\")\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Elaboration example\n",
    "\n",
    "- This part saves the spectra in fits files and updates some keywords"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from astropy.io import fits\n",
    "# This part saves the spectra in fits files and updates some keywords\n",
    "for year, data in data_by_ys.items():\n",
    "    print(year)\n",
    "    for ID,s in enumerate(data._p_list):\n",
    "        if (s.meta_data['src_name']==source_name):\n",
    "            if(s.meta_data['product']=='isgri_spectrum'):\n",
    "                ID_spec=ID\n",
    "            if(s.meta_data['product']=='isgri_arf'):\n",
    "                ID_arf=ID\n",
    "            if(s.meta_data['product']=='isgri_rmf'):\n",
    "                ID_rmf=ID\n",
    "\n",
    "    print(ID_spec, ID_arf, ID_rmf)\n",
    "\n",
    "    spec=data._p_list[ID_spec].data_unit[1].data\n",
    "    arf=data._p_list[ID_arf].data_unit[1].data\n",
    "    rmf=data._p_list[ID_rmf].data_unit[2].data\n",
    "    expos=data._p_list[0].data_unit[1].header['EXPOSURE']\n",
    "    name=source_name+'_'+str(year)\n",
    "    specname=name+'_spectrum.fits'\n",
    "    arfname=name+'_arf.fits.gz'\n",
    "    rmfname=name+'_rmf.fits.gz'\n",
    "    data._p_list[ID_spec].write_fits_file(specname)\n",
    "    data._p_list[ID_arf].write_fits_file(arfname)\n",
    "    data._p_list[ID_rmf].write_fits_file(rmfname)\n",
    "    hdul = fits.open(specname, mode='update')\n",
    "    hdul[1].header.set('EXPOSURE', expos)\n",
    "    hdul[1].header['RESPFILE']=rmfname\n",
    "    hdul[1].header['ANCRFILE']=arfname\n",
    "    hdul[1].data['SYS_ERR']=systematic_fraction\n",
    "\n",
    "    hdul.close()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Elaboration 2\n",
    "- If xspec is available, we make a fit of each spectrum"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "try:\n",
    "\n",
    "    import xspec\n",
    "    import shutil\n",
    "    from IPython.display import Image\n",
    "    from IPython.display import display \n",
    "\n",
    "    xspec.Fit.statMethod = \"chi\"\n",
    "\n",
    "    #init dictionaries\n",
    "    fit_by_lt={}\n",
    "\n",
    "    model='cflux*bknpow'\n",
    "\n",
    "    xspec.AllModels.systematic=0.0\n",
    "    low_energies=[20]\n",
    "    freeze_pow_ebreak=1\n",
    "\n",
    "    for year in range(start_year,end_year+1):\n",
    "\n",
    "        for c_emin in low_energies: #np.linspace(17,40,5):    \n",
    "            xspec.AllData.clear()\n",
    "\n",
    "            m1=xspec.Model(model)\n",
    "\n",
    "            specname=source_name+'_'+str(year)+'_spectrum.fits'\n",
    "\n",
    "            xspec.AllData(specname)\n",
    "\n",
    "            s = xspec.AllData(1)\n",
    "\n",
    "            isgri = xspec.AllModels(1)\n",
    "\n",
    "            print(m1.nParameters)\n",
    "\n",
    "            xspec.AllData.ignore('bad')\n",
    "            xspec.AllData.ignore('500.0-**')\n",
    "\n",
    "            ig=\"**-%.2f,500.-**\"%c_emin\n",
    "            print(\"ISGRI ignore: \"+ ig)\n",
    "            s.ignore(ig)\n",
    "\n",
    "            #Key for output\n",
    "            lt_key='%d_%.10lg'%(year, c_emin)\n",
    "\n",
    "            isgri.cflux.lg10Flux=-8            \n",
    "\n",
    "            isgri.cflux.Emin=20.\n",
    "            isgri.cflux.Emax=80.\n",
    "\n",
    "            isgri.bknpower.norm = \"1,-1\"\n",
    "            isgri.bknpower.PhoIndx1 = \"2.0,.01,1.,1.,3.,3.\"\n",
    "            isgri.bknpower.PhoIndx2 = \"2.2,.01,1.,1.,3.,3.\"\n",
    "            isgri.bknpower.BreakE = \"100,-1,20,20,300,300\"\n",
    "\n",
    "            xspec.Fit.perform()\n",
    "            isgri.bknpower.BreakE.frozen = freeze_pow_ebreak  > 0\n",
    "\n",
    "            xspec.Fit.perform()\n",
    "\n",
    "            max_chi=np.ceil(xspec.Fit.statistic / xspec.Fit.dof)\n",
    "\n",
    "            xspec.Fit.error(\"1.0 max %.1f 1-%d\"%(max_chi,m1.nParameters))\n",
    "\n",
    "\n",
    "            fit_by_lt[lt_key]=dict(\n",
    "                    emin=c_emin,\n",
    "                    year=year,\n",
    "                    chi2_red=xspec.Fit.statistic/xspec.Fit.dof,                                \n",
    "                    chi2=xspec.Fit.statistic,\n",
    "                    ndof=xspec.Fit.dof,                                    \n",
    "                )\n",
    "\n",
    "            for i in range(1,m1.nParameters+1): \n",
    "                if (not isgri(i).frozen) and (not bool(isgri(i).link)):\n",
    "                    #use the name plus position because there could be parameters with same name from multiple \n",
    "                    #model components (e.g., several gaussians)\n",
    "                    print(isgri(i).name, \"%.2f\"%(isgri(i).values[0]), isgri(i).frozen,bool(isgri(i).link) )\n",
    "                    fit_by_lt[lt_key][isgri(i).name+\"_%02d\"%(i)]=[ isgri(i).values[0], isgri(i).error[0], isgri(i).error[1] ]\n",
    "\n",
    "\n",
    "\n",
    "            xspec.Plot.device=\"/png\"\n",
    "            #xspec.Plot.addCommand(\"setplot en\")\n",
    "            xspec.Plot.xAxis=\"keV\"\n",
    "            xspec.Plot(\"ldata del\")\n",
    "            xspec.Plot.device=\"/png\"\n",
    "\n",
    "            fn=\"fit_%s.png\"%lt_key\n",
    "            fit_by_lt[lt_key]['plot_fname'] = fn\n",
    "\n",
    "            shutil.move(\"pgplot.png_2\", fn)\n",
    "\n",
    "            _=display(Image(filename=fn,format=\"png\"))\n",
    "\n",
    "except ImportError:\n",
    "    print(\"no problem!\")"
   ]
  }
 ],
 "metadata": {
  "celltoolbar": "Tags",
  "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.13.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}