{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Verify dust templates\n"
   ],
   "id": "f8cf04d96ecb"
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1a15c079-f15d-42fb-b66e-77843711d11a",
   "metadata": {},
   "outputs": [],
   "source": [
    "import os\n",
    "\n",
    "os.environ[\n",
    "    \"OMP_NUM_THREADS\"\n",
    "] = \"64\"  # for jupyter.nersc.gov otherwise the notebook only uses 2 cores"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9a40afd7-f8b9-45e1-a618-4a0f1053b925",
   "metadata": {},
   "outputs": [],
   "source": [
    "import pysm3\n",
    "import healpy as hp\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "12da229e-d783-4800-9378-4b9929d1e955",
   "metadata": {},
   "outputs": [],
   "source": [
    "cd /global/cfs/cdirs/cmb/www/pysm-data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c4602125-d810-41d5-9bf9-050cd33803ad",
   "metadata": {},
   "outputs": [],
   "source": [
    "cd dust_gnilc/"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5d0b279b-d58f-4767-8707-908141255b0c",
   "metadata": {},
   "outputs": [],
   "source": [
    "#name = \"template\"\n",
    "#name = \"Td\"\n",
    "name = \"beta\""
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9bce9499-e85e-4c35-afe9-9bc7db5191de",
   "metadata": {},
   "outputs": [],
   "source": [
    "pol = (0,1,2) if name == \"template\" else (0,)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "72bec66b-d27d-443d-af34-8a19ef70691c",
   "metadata": {},
   "outputs": [],
   "source": [
    "template_date = \"_2023.06.06\" if name in [\"beta\", \"Td\"] else \"\""
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "8c636ba1-7459-4f96-9089-1440e6e34987",
   "metadata": {},
   "outputs": [],
   "source": [
    "m = hp.read_map(f\"gnilc_dust_{name}_nside4096{template_date}.fits\", pol)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7d4f6a0c-a7d9-438e-87d9-a6548737e51f",
   "metadata": {},
   "outputs": [],
   "source": [
    "m2048 = hp.read_map(f\"gnilc_dust_{name}_nside2048{template_date}.fits\", pol)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7eb53f5e-a7d4-4090-8f13-08faee35bd61",
   "metadata": {},
   "outputs": [],
   "source": [
    "m8192 = hp.read_map(f\"gnilc_dust_{name}_nside8192{template_date}.fits\", pol)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "953817f1-39bf-49f9-8947-b3d503514d42",
   "metadata": {},
   "outputs": [],
   "source": [
    "cl = hp.anafast(m, lmax=int(2.5*4096)).reshape((len(pol), -1))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "eab9af0e-9cb7-430c-b4b5-4399dc58f1b9",
   "metadata": {},
   "outputs": [],
   "source": [
    "cl8192 = hp.anafast(m8192, lmax=int(2*8192)).reshape((len(pol), -1))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f8af6555-0863-4db2-9e92-9c963e945d0b",
   "metadata": {},
   "outputs": [],
   "source": [
    "cl2048 = hp.anafast(m2048, lmax=int(2.5*2048)).reshape((len(pol), -1))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "b4893273-385c-4ef0-aa13-4feaedd159f9",
   "metadata": {},
   "outputs": [],
   "source": [
    "from pathlib import Path"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "36d90f03-08f7-4dc2-adac-0275d2156858",
   "metadata": {},
   "outputs": [],
   "source": [
    "datadir = Path(\"raw\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "971ec2f6-94dc-4f0b-ab6a-cceb80f9f1ec",
   "metadata": {},
   "outputs": [],
   "source": [
    "if name == \"template\":\n",
    "    name = \"logpoltens\"\n",
    "    suffix=\"_complex64\"\n",
    "else:\n",
    "    suffix = \"\""
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "b5fdb5bf-7cb1-4751-a634-c6027ebdad07",
   "metadata": {},
   "outputs": [],
   "source": [
    "alm_large_scale = hp.read_alm(\n",
    "    datadir / f\"gnilc_dust_largescale_template_{name}_alm_nside2048_lmax1024{suffix}.fits.gz\",\n",
    "    hdu=1,\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "2003da95-0d7b-49bc-a96e-ec370dacc1bd",
   "metadata": {},
   "outputs": [],
   "source": [
    "cl_small_scale = hp.read_cl(\n",
    "    datadir / f\"gnilc_dust_small_scales_{name}_cl_lmax16384{template_date}.fits.gz\"\n",
    ").reshape((len(pol), -1))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "cf32be9b-f75f-4b14-8ab9-4d7a0aa94c99",
   "metadata": {},
   "outputs": [],
   "source": [
    "if name == \"logpoltens\":\n",
    "    name = \"template\""
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "61b405a9-6aeb-415b-9657-8aca19e37007",
   "metadata": {},
   "outputs": [],
   "source": [
    "pol_label = \"TEB\""
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6cb6aae8-8fa8-460d-a0ab-9915b26e6aae",
   "metadata": {},
   "outputs": [],
   "source": [
    "for p in pol:\n",
    "    plt.figure(figsize=(12,6))\n",
    "    plt.loglog(cl8192[p], label=f\"map at Nside=8192\")\n",
    "    plt.loglog(cl[p], label=\"map at Nside=4096\")\n",
    "    plt.loglog(cl2048[p], label=f\"map at Nside=2048\")\n",
    "    if name != \"template\":\n",
    "        plt.loglog(hp.alm2cl(alm_large_scale.astype(complex).reshape((len(pol), -1)))[p], \"--\", label=\"large scale\")\n",
    "        plt.loglog(cl_small_scale[p], \"--\", label=\"small scales\")\n",
    "    plt.legend()\n",
    "    plt.title(f\"{name} maps spectra comparison {pol_label[p]}\")\n",
    "    plt.grid();"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "cfafe7e1-6445-4611-8678-dcb33c7a62f4",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "language_info": {
   "codemirror_mode": {
    "name": "ipython"
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}