{ "cells": [ { "cell_type": "markdown", "id": "6fa69231-aae2-4cb6-b131-30864c5f2e3a", "metadata": {}, "source": [ "# E: Ejercicio resuelto: Generacion de una funcion para reducir datos" ] }, { "cell_type": "markdown", "id": "6310b5b4-0fda-4638-864d-0642a5938698", "metadata": {}, "source": [ "***\n", "**Ejercicio E.1:**\n", "\n", "Basándonos en el código de los anteriores `notebooks` crear una función a la que pasándole el directorio donde están nuestras imágenes, bias, darks y flats reduzca la imagen automáticamente." ] }, { "cell_type": "code", "execution_count": 1, "id": "5a709831-eb65-4934-af5d-8bb0a374c7a2", "metadata": {}, "outputs": [], "source": [ "import matplotlib.pyplot as plt\n", "\n", "import numpy as np\n", "import glob\n", "from astropy.io import fits\n", "import os\n", "import astroscrappy" ] }, { "cell_type": "code", "execution_count": 2, "id": "4953a72b-8ca6-4dff-8816-9ada9a2ff092", "metadata": {}, "outputs": [], "source": [ "import astroalign\n", "from astropy.io import fits as ft\n", "\n", "\n", "def imprimir_info(p, ii):\n", " # Esta función imprime por pantalla la info de la transformación que se aplicará.\n", " # Aunque es innecesaria, sirve para que el usuario sepa que la máquina está haciendo algo.\n", " print(\"\\nAlineando imagen {:}\".format(ii))\n", " print(\"Rotación: {:.2f} grados\".format(p.rotation * 180.0 / np.pi))\n", " print(\"Factor de escala: {:.2f}\".format(p.scale))\n", " print(\"Traslación: (x, y) = ({:.2f}, {:.2f})\".format(*p.translation))\n", "\n", "\n", "def registra_lista(lista, infoAlineado = False):\n", " '''\n", " Función para la alineacion de imagenes a traves de astroalign\n", " \n", " lista: tenemos que pasar una lista de la ruta donde se encuentran nuestras imagenes. List of string.\n", " \n", " Return: Lista de numpy array con los datos de las imagenes alineadas. List of numpy array. \n", " '''\n", " \n", " \n", " cantidad=len(lista)\n", " \n", " # La primera imagen de la lista será la toma de referencia.\n", " \n", " print(\"\\nComenzando la alineación.\")\n", " print(\"\\nLa toma de referencia es {:}\".format(lista[0])) \n", " blanco=ft.open(lista[0])\n", " img_blanco=blanco[0].data\n", " hdr_blanco=blanco[0].header\n", " blanco.close()\n", "\n", " del(lista[0]) # Quito la imagen de referencia del listado\n", " imagenes_alineadas = [img_blanco, ]\n", " for ii in lista:\n", " ff=ft.open(ii)\n", " img_torcida=ff[0].data\n", " hdr_torcida=ff[0].header\n", " ff.close()\n", " p, (pos_img, pos_img_rot) = astroalign.find_transform(img_torcida, img_blanco)\n", " img_aligned = astroalign.register(img_torcida, img_blanco)\n", " hdr_torcida.add_comment(\"Registrado con Astroalign y PyReduc\")\n", " imagenes_alineadas.append(img_aligned[0])\n", " if infoAlineado:\n", " imprimir_info(p, ii)\n", "\n", " print(\"\\nAlineado realizado con éxito\")\n", " return imagenes_alineadas" ] }, { "cell_type": "markdown", "id": "e5a83d5f-692a-4779-912f-e2d108a5759f", "metadata": {}, "source": [ "Creamos nuestra funcion reducir en donde las entradas de la función serán los directorios donde tenemos guardados los archivos `FITS`." ] }, { "cell_type": "code", "execution_count": 3, "id": "30250cc3-1c57-4799-a455-1f2295393a24", "metadata": {}, "outputs": [], "source": [ "def reducir(dirBias, dirDarks, dirFlats, dirLights, dirSalida, imaInter = False, infoCosmicray = False):\n", " '''\n", " Función para reducir datos en formato fit. En los directorios solo deben\n", " estar los ficheros particulares de cada caso.\n", " \n", " dirBias: directorio donde se encuentran los bias. String.\n", " dirDarks: directorio donde se encuentran los darks. String.\n", " dirFlats: directorio donde se encuentran los bias. String.\n", " dirLights: directorio donde se encuentran los bias. String.\n", " dirSalida: directorio donde irán las imagenes calibradas en formato fit. String.\n", " imaInter: Si queremos guardar las imagenes intermedias pondremos True. Booleano.\n", " infoCosmicray: Si queremos informacion de los rayos cosmicos que ha detectado el algoritmo astroscrappy. Booleano.\n", " \n", " Return: Imagen calibrada en forma de matriz. numpy array. \n", " '''\n", "\n", " try:\n", " os.mkdir(dirSalida)\n", " except:\n", " print(f\"El directorio {dirSalida} ya existe\")\n", " \n", " bias_list = sorted(glob.glob(f\"{dirBias}*\"))\n", " darks_list = sorted(glob.glob(f\"{dirDarks}*\"))\n", " flats_list = sorted(glob.glob(f\"{dirFlats}*\"))\n", " lights_list = sorted(glob.glob(f\"{dirLights}*\"))\n", " \n", " # Bias\n", " \n", " print(\"Procesando BIAS\")\n", " \n", " bias_data = []\n", " bias_header = []\n", " for bias in bias_list:\n", " hdul_bias = fits.open(bias)\n", " bias_data.append(hdul_bias[0].data) \n", " bias_header.append(hdul_bias[0].header)\n", " hdul_bias.close()\n", " for r in range(len(bias_list)):\n", " if np.nanstd(bias_data[r]) > np.median(np.nanstd(bias_data))*1.2:\n", " print(\"Eliminando BIAS: \" + bias_list[r])\n", " bias_list.remove(bias_list[r])\n", " superBias = np.nanmedian(bias_data, axis=0)\n", " bias_header[0]['history'] = f\"Esta imagen es un stacking de {len(bias_list)} imagenes\" \n", " \n", " # Darks\n", " \n", " print(\"Procesando DARKs\")\n", " \n", " darks_data = []\n", " darks_header = []\n", " darks_current = []\n", " for dark in darks_list:\n", " hdul_darks = fits.open(dark)\n", " darks_data.append(hdul_darks[0].data)\n", " darks_header.append(hdul_darks[0].header)\n", " hdul_darks.close()\n", " for dark in range(len(darks_list)):\n", " darks_current.append(darks_data[dark] - superBias)\n", " for i in range(len(darks_list)):\n", " darks_current[i][darks_current[i] < 0] = 0\n", " superDark = np.nanmedian(darks_data, axis=0)\n", " darkCurrent = np.median(darks_current, axis=0)\n", " \n", " #Flats\n", " \n", " print(f\"Procesando FLATs\")\n", " \n", " flats_data = []\n", " flats_header = []\n", " flat_data_b = []\n", " for flat in flats_list:\n", " hdul_flats = fits.open(flat)\n", " flats_data.append(hdul_flats[0].data) \n", " flats_header.append(hdul_flats[0].header)\n", " hdul_flats.close()\n", " for r in range(len(flats_list)):\n", " if np.nanstd(flats_data[r]) > np.median(np.nanstd(flats_data))*1.1:\n", " print(bias_list[r])\n", " bias_list.remove(bias_list[r])\n", " y_len, x_len = np.shape(flats_data[0])\n", " x_center, y_center = x_len//2, y_len//2\n", " for flat in range(len(flats_list)):\n", " flat_superBias = flats_data[flat] - superBias\n", " flat_norm = flat_superBias / np.median(flat_superBias[y_center-250:y_center+250,x_center-250:x_center+250])\n", " flat_data_b.append(flat_norm)\n", " superFlat = np.nanmedian(flat_data_b, axis=0)\n", " \n", " #lights\n", " \n", " print(f\"Procesando LIGHTs de {flats_header[0]['FILTER']}\")\n", " \n", " lights_header = []\n", " lights_data = []\n", " FWHM = []\n", " for f in lights_list:\n", " hdul = fits.open(f)\n", " lights_header.append(hdul[0].header)\n", " lights_data.append(hdul[0].data)\n", " FWHM.append(hdul[0].header['FWHM'])\n", " hdul.close()\n", "# print(f\"Lights de {lights_header[0]['FILTER']}\")\n", " for r in range(len(lights_list)):\n", " if lights_header[r]['FWHM'] > np.median(FWHM)*1.5:\n", " lights_data.remove(lights_data[r])\n", " lights_header.reomve(lights_header[r])\n", " lights_list.remove(lights_list[r])\n", " print(f\"Se ha eliminado la imagen {lights_list[r]}\")\n", " \n", " lights_correction = []\n", " for light in range(len(lights_list)):\n", " light_b =lights_data[light].astype(np.float32) - superBias.astype(np.float32)\n", " light_bd = light_b - darkCurrent.astype(np.float32) \n", " light_bdf = light_bd / superFlat.astype(np.float32) \n", " lights_correction.append(light_bdf)\n", " \n", " for hot in range(len(lights_list)):\n", " hot_pixeles = np.where(lights_correction[hot] < 0) # hotpixel será un array con la posicion del los hot pixeles\n", " for i in range(len(hot_pixeles[0])): # \n", " lights_correction[hot][hot_pixeles[0][i]][hot_pixeles[1][i]] = np.median([lights_correction[hot][hot_pixeles[0][i]-1, hot_pixeles[1][0]], lights_correction[hot][hot_pixeles[1][i]+1, hot_pixeles[1][0]]]) \n", " lights_correction_cosmicray = []\n", " mask_cosmicray = []\n", " for image in range(len(lights_list)):\n", " mask, cosmicray = astroscrappy.detect_cosmics(lights_correction[image], sigclip=6, sigfrac=15, objlim=5.0, gain=1.0, readnoise=6.5, satlevel=65536.0, niter=4, verbose=infoCosmicray)\n", " lights_correction_cosmicray.append(cosmicray)\n", " mask_cosmicray.append(mask)\n", " imagenes_alineadas = registra_lista(lights_list)\n", " imagenFinal = np.nanmedian(imagenes_alineadas, axis=0)\n", " # Guardamos\n", " if imaInter:\n", " hdu_Bias = fits.PrimaryHDU(data=superBias.astype(np.float32), header=bias_header[0])\n", " hdu_Bias.writeto(dirSalida + \"superBias.fit\", overwrite=True)\n", "\n", " super_dark = fits.PrimaryHDU(data=superDark.astype(np.float32), header=darks_header[0])\n", " super_dark.writeto(dirSalida + \"superDark.fit\", overwrite=True)\n", " hdu_Dark_current = fits.PrimaryHDU(data=darkCurrent.astype(np.float32))\n", " hdu_Dark_current.writeto(dirSalida + \"darkCurrent.fit\", overwrite=True)\n", "\n", " hdu_FlatH = fits.PrimaryHDU(data=superFlat.astype(np.float32), header=flats_header[0])\n", " hdu_FlatH.writeto(dirSalida + \"superFlat.fit\",overwrite=True)\n", " \n", " hdu_imagenFinal = fits.PrimaryHDU(data=imagenFinal.astype(np.float32), header=lights_header[0])\n", " hdu_imagenFinal.writeto(dirSalida + f\"imagenFinal_{lights_header[0]['FILTER']}.fit\",overwrite=True)\n", " return imagenFinal" ] }, { "cell_type": "markdown", "id": "ff7561de-9a79-4aa6-961b-13522e72632b", "metadata": {}, "source": [ "Generamos los directorios y llamamos a la funcion ``reducir``." ] }, { "cell_type": "code", "execution_count": 4, "id": "6878bad9-e082-400e-b344-729e29e65034", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Probablemente el directorio ya existe\n", "El directorio salidas/salidaCalibracion/sho/ ya existe\n", "Procesando BIAS\n", "Eliminando BIAS: imagenes/calibracionImagenes/bias/Telescope3_BIAS_1x1_0s_-30degC_20200120_0547_000001089.fit\n", "Procesando DARKs\n", "Procesando FLATs\n", "Procesando LIGHTs de H-alpha\n", "\n", "Comenzando la alineación.\n", "\n", "La toma de referencia es imagenes/calibracionImagenes/ha/NGC2070-20200208@030339-300S-H-alpha.fit\n", "\n", "Alineado realizado con éxito\n", "El directorio salidas/salidaCalibracion/sho/ ya existe\n", "Procesando BIAS\n", "Eliminando BIAS: imagenes/calibracionImagenes/bias/Telescope3_BIAS_1x1_0s_-30degC_20200120_0547_000001089.fit\n", "Procesando DARKs\n", "Procesando FLATs\n", "Procesando LIGHTs de Oiii\n", "\n", "Comenzando la alineación.\n", "\n", "La toma de referencia es imagenes/calibracionImagenes/O3/NGC2070-20200208@035232-300S-Oiii.fit\n", "\n", "Alineado realizado con éxito\n", "El directorio salidas/salidaCalibracion/sho/ ya existe\n", "Procesando BIAS\n", "Eliminando BIAS: imagenes/calibracionImagenes/bias/Telescope3_BIAS_1x1_0s_-30degC_20200120_0547_000001089.fit\n", "Procesando DARKs\n", "Procesando FLATs\n", "Procesando LIGHTs de Sii\n", "\n", "Comenzando la alineación.\n", "\n", "La toma de referencia es imagenes/calibracionImagenes/S2/NGC2070-20200209@020725-600S-Sii.fit\n", "\n", "Alineado realizado con éxito\n" ] } ], "source": [ "import os\n", "\n", "baseDir = \"imagenes/calibracionImagenes/\"\n", "\n", "dirBias = baseDir + 'bias/' \n", "dirDarks = baseDir + 'darks/'\n", "\n", "dirFlatsHa = baseDir + 'flats/ha/'\n", "dirFlatsO3 = baseDir + 'flats/O3/'\n", "dirFlatsS2 = baseDir + 'flats/S2/'\n", "\n", "dirLightsHa = baseDir + 'ha/' \n", "dirLightsO3 = baseDir + 'O3/'\n", "dirLightsS2 = baseDir + 'S2/'\n", "\n", "dirSalida = 'salidas/salidaCalibracion/sho/'\n", "\n", "try:\n", " os.mkdir(dirSalida)\n", "except:\n", " print(\"Probablemente el directorio ya existe\")\n", "\n", "ha = reducir(dirBias, dirDarks, dirFlatsHa, dirLightsHa, dirSalida, imaInter = False)\n", "O3 = reducir(dirBias, dirDarks, dirFlatsO3, dirLightsO3, dirSalida, imaInter = False)\n", "S2 = reducir(dirBias, dirDarks, dirFlatsS2, dirLightsS2, dirSalida, imaInter = False)" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "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.16" } }, "nbformat": 4, "nbformat_minor": 5 }