REDUCCION DE IMÁGENES CCD CON IRAF


      Descripción de pasos en la reducción
El proceso de reducción implica la eliminación de efectos instrumentales que estén presentes en los datos, ya sean de espectroscopía o de imagen directa, y que es necesario llevarlo a cabo antes de poder realizar cualquier tipo de medida sobre nuestros datos.

Solamente nos ocuparemos de los efectos más comunes y que aparecen en todos los datos CCD:

    Determinación del nivel de pedestal o nivel cero ("bias") - Se trata de eliminar un nivel de cero que contienen los datos y que es añadido electrónicamente, simplemente para mantener la positividad del numero medido. Este valor puede variar pixel a pixel a través del detector y también a lo largo de la noche.
    Para determinar la estructura espacial se toman series de exposiciones con tiempo de integración cero al principio o final de la noche.
     
      ¿Cuántas imágenes de cero se deben obtener?

      En primer lugar debemos decidir si queremos sustraer una imagen de cero, para ello es necesario estimar cual es el tipo de ruido en nuestros datos y hacer un compromiso entre el nivel de variación observado en la imagen de cero y el nivel de ruido de la imagen resultante. Si estamos en el limite donde domina el ruido de lectura en las imágenes del cielo (por ejemplo exposiciones cortas, filtros de banda estrecha, o espectroscopia de alta resolución), tendremos que evitar añadir más ruido al sustraer una imagen de cero. La imagen de nivel cero se genera a partir de un promedio de N_zero imagenes tomadas con tiempo de integración cero. El ruido de la imagen resultante, despues de sustraer la imagen cero promedio, se puede calcular aproximadamente por la siguiente expresión:
       

          Ruido_final = Rd * SQRT(1+1/N_zero),


      siendo Rd el ruido de lectura del CCD.
      Así pues si tenemos sólo una imagen de cero para promediar aumentaremos el ruido un factor SQRT(2), mientras que si promediamos 10 imágenes de cero añadiremos solo un 5% más de ruido. Este valor es el que debemos comparar con las variaciones  observadas sobre la imagen de cero.


    El nivel de cero puede cambiar a lo largo de la noche, por lo cada imagen suele tener una zona donde se almacena el nivel de cero propio, la llamada zona de "overscan". Esta zona está formada por unas columnas o filas adicionales, que suelen contener todas las imágenes. Estas filas o columnas son virtuales, no existen físicamente sobre el detector, siendo creadas artificialmente por la electrónica de lectura, y por tanto no están iluminadas. Sin embargo, este nivel no es representativo de todos los pixeles, sino únicamente de unas pocas columnas o filas. Por tanto se suele
    promediar y obtener solo un valor por imagen o en algunos casos un ajuste lineal a lo largo de las columnas o filas de "overscan".

     
    Corrección de sensibilidad ("flat-field'') - Se trata de corregir la no uniformidad de la respuesta pixel a pixel y a gran escala. Para ello se necesitan imágenes de un campo suficientemente plano y con un nivel alto de cuentas para reducir el ruido (>10000). En general se pueden tener tres tipos de "flat-fields": del propio cielo, que se toman al atardecer o amanecer; de cúpula iluminada o iluminación directa por lampara halógena. Los primeros son sin duda los más recomendados pero también son más difíciles de obtener, sobre todo cuando se hace espectroscopía o imagen en filtros de banda estrecha por la rápida variación del brillo del cielo.
     
      ¿Cuantas imágenes de "flat'' se deben obtener?

      La corrección de flat-field debe servir como minimo para eliminar las variaciones en sensibilidad de pixel a pixel, por ello es necesario obtener exposiciones de un campo relativamente uniforme. La precisión de estas correcciones influye de forma definitiva sobre la precisión de las medidas que se hagan a posteriori sobre la imagen, ya sea fotometría o espectroscopía. Así pues si queremos tener un error del orden del 1% comparando pixel a pixel, es necesario que nuestra imagen de flat-field no supere ese valor. Para ello es necesario tomar imágenes cuyo valor promedio de cuentas sea alto. Hay que recordar que la relación señal/ruido en el limite que estamos, dominado por ruido fotónico y para el promedio de varias imágenes, está dada por la siguiente expresión:
       

        S/N = SQRT(S_ADU*G*N_aver)


      Siendo G el factor de conversion de electrones a ADUs (usualmente en el rango 1 -2), S_ADU el valor promedio de la imagen, y N_aver el numero de imágenes promediadas. Así pues para mejorar la relación S/N y tener un error por debajo del 1% es necesario que S_ADU sea mayor que 10000. Por otra parte hay que tener en cuenta no saturar el detector y no tomar valores excesivamente altos dado que tendremos problemas de no linealidad. En general se toman varias imágenes para limpiar posibles impactos de rayos cósmicos, siendo 5 o más un número adecuado para este propósito.
       
       

    Pixeles muertos o no lineales - Dependiendo de la calidad del CCD se pueden tener más o menos "pixeles malos". Estos pixeles son aquellos que no responden como la mayoria a la luz que reciben o más comúnmente no tienen una respuesta lineal. En algunas aplicaciones se suelen corregir por interpolación de los pixeles vecinos (por ejemplo de forma interactiva usando imedit).
       
    Determinación de sección overscan y sección útil del CCD

    El primer paso antes de reducir los datos es determinar cual es la sección de la imagen donde se encuentra la zona de "overscan", a la cual nos referiremos como "biassec". Además hay que determinar la sección que contiene datos buenos, es decir sin viñeteo y fuera de la seccion de "overscan", a esta sección se la denomina "trimsec".

    Para esto se examina una buena imagen de "flat-field" con el programa implot. Haremos un promedio (":a 15", promedia usando 15 filas o columnas) por filas (":l 520", dibuja la fila 520) y por columnas (":c 600", dibuja la columna 600) para ver hasta donde hay una iluminación uniforme por los extremos. Así se determina la zona con datos útiles, o "trimsec". La zona de "overscan" se encuentra fuera de la zona iluminada del CCD, y se reconoce porque tiene el nivel más bajo de cuentas y no muestra estructura como la parte iluminada. A esta sección de la imagen se denomina "biassec".

    Recientemente algunos instrumentos colocan en las cabeceras de los datos que producen, las dos palabras clave biassec y trimsec. Como buena práctica se recomienda comprobar que esa información es correcta.

    Hay que tener en cuenta que IRAF usa la notación sec_imagen[x1:x2,y1:y2], para describir la parte de la imagen que va desde la columna x1 a la columna x2, y desde la fila y1 a la fila y2.
     

       
      Sustracción de nivel cero

      Para sustraer el nivel de cero impuesto por la electrónica de las imagenes ya hemos visto que existen dos posibilidades:

        Si se han tomado imágenes con tiempo de integración cero, se pueden combinar estas mismas y obtener una imagen que representa el nivel cero promedio para cada pixel.
        El nivel de cero puede cambiar a lo largo de la noche, por lo cada imagen suele tener una zona donde se almacena el nivel de cero propio, la llamada zona de "overscan". Sin embargo, este nivel no es representativo de todos los pixeles, sino únicamente de unas pocas columnas o filas. Por tanto se suele
        promediar y obtener solo un valor por imagen o en algunos casos un ajuste lineal a lo largo de las columnas o filas de "overscan".
         
        Combinando imágenes de zero

        En caso que tenemos suficientes imágenes de cero se pueden combinar con el programa zerocombine. A continuación se muestran sus parámetros de defecto:


         

        Notar que se obtiene el valor promedio (combine="average"), descartando el valor máximo para cada pixel (reject="minmax"; nlow=0; nhigh=1).
         

        Primer paso por ccdproc


      El programa ccdproc es el encargado de realizar la mayor parte del tratamiento de las imágenes CCD, desde sustraer el nivel de cero, hasta corregir flat-field y pixeles malos. El primer paso a través de ccdproc sirve para cortar las imagenes, sustraer un nivel de cero constante para cada imagen, y en su caso, sustraer también la imagen Zero que contiene la estructura espacial del nivel de cero. A continuación se muestra la lista de los parámetros de ccdproc para hacer este procesamiento:


       
       

      Corrección de respuesta pixel a pixel ("flatfielding")

         

        Combinando imágenes de flat-field

        Una vez hemos visto que imágenes de flat-field tienen un nivel de cuentas adecuado, se trata de combinarlas para obtener una imagen limpia. Hay que tener en cuenta que se combinaran en el caso de imagen directa solo aquellas que corresponden al mismo filtro y en el caso de espectroscopia aquellas en las que el rango espectral cubierto sea aproximadamente el mismo. Para combinar imágenes de flat-field se usa el programa flatcombine. A continuación se muestra un ejemplo de selección de parámetros:

        En este caso se obtiene el resultado (combine="average") como el promedio de todos los valores, antes de obtener el valor final promedio se eliminan aquellos valores que exceden por un factor la desviacion tipica (reject="sigclip"; lsigma=3; hsigma=3). Además antes de combinar las distintas imágenes se escalan por el valor de la moda (scale="mode") de la distribución de todos sus pixeles para que todas tengan un valor medio de cuentas similar. En el ejemplo anterior, el resultado de flatcombine será un conjunto de imágenes de flat-field, una por cada uno de los filtros que se hayan observado. Para que esto funcione correctamente es necesario que el parametro flatcombine.subsets=yes, además imagetyp y subset deben estar definidos correctamente a flat y la palabra clave en la cabecera que distingue el filtro). Más simplemente se puede usar flatcombine separadamente para cada conjunto de imágenes, para ello se formaran listas con las imágenes de flat-field que corresponden al mismo filtro o rango en longitud de onda.
         

        Caso espectroscopía:

        En el caso de espectroscopia las imágenes de flat-field presentan estructuras muy diferentes a lo largo de la dirección espacial y la espectral. En general usaremos exposiciones de una lampara de luz blanca o la cúpula iluminada. A lo largo de la dirección espectral la estructura observada es la convolución de varias funciones: el color de la luz usada, la transmisión de la red de dispersión, y la respuesta espectral del detector. Lo que nos interesa corregir son las variaciones pixel a pixel que se presentan para cualquier tipo de iluminación, por tanto debemos eliminar esa estructura a gran escala que es peculiar al tipo de exposición y que no tiene porque aparecer en los objetos de nuestro interés. Para eliminar esta estructura se ajusta una función suave al promedio de todas las filas o columnas a través de la dirección espacial, y finalmente se determina la corrección a partir del cociente entre las dos. El programa que realiza esta tarea se llama response y se encuentra dentro del paquete twodspec y a su vez dentro de longslit. Por ejemplo:

        response flat flat nflat

        siendo los dos primeros argumentos, la imagen de calibración y la de normalización, flat (combinación de todos los flat-field) y el tercer argumento es el resultado, que será una imagen normalizada.

        Una vez dentro del programa response, estaremos usando uno de los entornos más usados dentro de IRAF, el paquete icfit, que se usa para realizar ajustes de funciones a datos en una dimensión. Los comandos que se usan comúnmente aquí son:
         

          :function [legendre|chebyshev|spline3|spline1]

          :order n- el orden es el número de términos, en polinomios es el grado más 1

          f - sirve para hacer un nuevo ajuste

          k - sirve para ver el cociente entre el espectro y su ajuste

          q - salir del programa y producir la imagen de flat-field normalizado.


         

        Normalmente es necesario usar funciones tipo spline3 y un numero de términos elevado para ajustar adecuadamente la variación espectral observada. Se debe procurar eliminar toda la variación observada a gran escala en la dirección espectral.
         

        Segundo paso por ccdproc


      Una vez producidas las imágenes de flat-field para cada filtro o rango espectral ya se pueden corregir las imágenes del programa de observación. Para ello crearemos una lista con todas ellas, separandolas por filtro si es necesario. El programa que aplica la corrección es de nuevo ccdproc, por ejemplo:
       

        ccdproc @lista_filtroV flatcor+ flat="FlatV"


      Hay que notar que en el caso de imagen directa no hemos normalizado la imagen de flat-field, sin embargo ccdproc se encargará de ello.
      Para comprobar que efecto ha tenido la corrección de flat-field, seleccionamos una imagen con tiempo de exposición relativamente largo para que el ruido que domine sea de tipo Poisson. Luego comparamos el ruido, o desviación estandar, en las dos imágenes (sin corregir y corregida) en una zona donde no haya objetos, se debe notar una disminución del mismo en la imagen que ha sido corregida, lo cual significa que hemos igualado la respuesta de los pixeles, lo cual mejorará la calidad de las medidas posteriores.
       

      Corrección de iluminación

      Una vez realizada la corrección de flat-field pixel a pixel pasaremos a comprobar que no quedan gradientes de iluminación a gran escala. Especialmente en el caso de espectroscopia que hemos usado imágenes de flat-field obtenidas con iluminación diferente a los objetos sobre el cielo (exposiciones de lampara o cúpula). Usualmente se comprueban exposiciones del cielo brillante (al atardecer o amanecer) ya corregidas por flat-field, promediando varias filas o columnas a través de la dirección espectral y buscando gradientes (implot es una buena herramienta para esto).

      En caso que sea necesario corregir este efecto, usaremos el programa illumination para crear una imagen correctora a partir de una imagen de cielo con suficientes cuentas (o combinación de ellas), ya corregida por el flat-field mencionado anteriormente. Este programa ajusta una función a lo largo de la dirección espacial sobre un promedio de filas obtenidas promediando a través de la dirección espectral, en varios rangos de longitud de onda.

      En el caso de imagen directa se puede crear esta misma corrección si fuese necesaria usando el programa mkskycor para suavizar imágenes de cielo con tiempos de exposición largos, eliminando estructuras a pequeña escala como estrellas.
       
       

        Tercer paso por ccdproc


      Para aplicar la correccion de iluminación a todas las exposiciones de nuestros objetos se usa de nuevo ccdproc, por ejemplo:
       

        ccdproc @lista_objetos illumcor+ illum="skyflat"
      Corrección de pixeles defectuosos
      Una forma muy eficiente para crear una imagen de pixeles malos es a partir de dos conjuntos de imágenes de flat-field con dos niveles de cuentas, uno alto y otro más bajo. Después de combinarlas por separado se divide una por otra y los pixeles malos aparecen claramente como desviaciones sobre la media. El programa ccdmask (dentro del paquete nmisc) construye una imagen máscara con valor 0 para los pixeles buenos, 1 o 2 para filas y columnas malas, respectivamente.

      El programa fixpix (dentro del paquete nmisc) sirve para corregir los pixeles malos registrados en la imagen mascara, por interpolación sobre filas o columnas dependiendo del valor en la mascara.


Preparación para una reducción comoda. Pasos previos a la automatización: setupinstrument, ccdlist
 

La reducción de datos CCD dentro de IRAF se realiza básicamente con el paquete ccdred. La mayoría de programas dentro de este paquete basan su funcionamiento en la información que encuentran en la cabecera de las imágenes. Por tanto la primera operación que se debe realizar es comprobar que la información contenida en las cabeceras es interpretada correctamente por los programas, y estos distinguen las imágenes de "bias", de aquellas de "flat-field", y de los propios objetos.
 
 

    Setupinstrument

    Necesitamos pues crear una correspondencia entre las palabras clave escritas en la cabecera y sus valores a aquellos que entiende el paquete ccdred. Para esto se usa el comando setinstrument, hay que cargar ccdred y imred. Se ejecutará de la forma siguiente:

    setinstrument direct site='iac80' directory='miccddb$' rev+

    Se asume que hemos definido la variable "miccddb" a un directorio que contiene a su vez un subdirectorio llamado "iac80". En ese directorio se creara un fichero direct.cl, que permite predefinir algunos parámetros como BIASSEC y TRIMSEC:

    # Generic routine for setting parameters.

    ccdred.pixeltype = "real real"

    ccdred.verbose = yes

    ccdred.logfile = "logfile"

    ccdred.plotfile = ""

    ccdred.backup = ""

    ccdred.instrument = "miccddbloc$iac80/direct.dat"

    ccdred.ssfile = "subsets"

    ccdred.graphics = "stdgraph"

    ccdred.cursor = ""

    ccdproc.fixpix = yes

    ccdproc.overscan = yes

    ccdproc.trim = yes

    ccdproc.zerocor = no

    ccdproc.darkcor = no

    ccdproc.flatcor = no

    ccdproc.readcor = no

    ccdproc.scancor = no

    ccdproc.readaxis = "line"

    ccdproc.biassec = "[1046:1055,2:1024]"

    ccdproc.trimsec = "[3:1002,25:1024]"

    ccdproc.interactive = no

    ccdproc.function = "legendre"

    ccdproc.order = 1

    ccdproc.sample = "*"

    ccdproc.naverage = 1

    ccdproc.niterate = 1

    ccdproc.low_reject = 3

    ccdproc.high_reject = 3

    ccdproc.grow = 0

    ccdproc.fixfile='miccddbloc$iac80/badpix.dat'
     
     

    También es necesario crear un fichero llamado "direct.dat" que contendrá la conversión entre "keywords", a continuación se puede ver un ejemplo adecuado al IAC80:

    subset insfilte

    DARK dark

    BIAS zero

    OBJECT object

    'DOME FLAT' flat

    'PROJECTOR FLAT' flat

    'SKY FLAT' object

    La primera columna corresponde a las palabras que espera ccdred, mientras que la columna de la derecha corresponde a lo que se encuentra en la cabecera. La palabra clave seleccionada como subset (en este caso insfilte) permite que los programas que combinan muchas imágenes lo hagan teniendo en cuenta el valor de la misma, que en este caso corresponde al filtro seleccionado.
     

    Agrupando imágenes según tipos de imágenes

Una vez hemos creado los ficheros de traslación y ejecutado el comando setinstrument, se puede comprobar que IRAF reconoce correctamente los distintos tipos de imágenes (bias, flat-field, objects, etc). Para esto se usa el comando ccdlist sobre todas las imágenes. Obtendremos un resultado similar al siguiente:
cc>ccdlist feb*.fits ccdtype='object'
feb150049.fits[1056,1024][ushort][object][]:Rubin149_B 120s
feb150050.fits[1056,1024][ushort][object][]:Rubin149_V 120s
feb150051.fits[1056,1024][ushort][object][]:Rubin149_R 90s
feb150061.fits[1056,1024][ushort][object][]:NGC2713_R
feb150062.fits[1056,1024][ushort][object][]:NGC2713_B
feb150068.fits[1056,1024][ushort][object][]:Rubin149_V
feb150069.fits[1056,1024][ushort][object][]:Rubin149_R
feb150070.fits[1056,1024][ushort][object][]:Rubin149_B
feb150082.fits[1056,1024][ushort][object][]:SAO104430_R
feb150083.fits[1056,1024][ushort][object][]:SAO104430_V
feb150084.fits[1056,1024][ushort][object][]:SAO104430_B
      Referencias