The SIM class (structured images)
in the astronomy & astrophysics toolbox for MATLAB
Description:
The SIM class is designed to store, access and manipulate images and their metadata.
This class contains a large number of methods to analyze and reduce astronomical images, from the most basic operations (e.g., arithmatic, bias, flat) to high-level functions (e.g., astrometry, registration and image subtraction).
The SIM class inherits from the WorldCooSys, HEAD, AstCat, MASK, and ClassPSF classes.
This file is accessible through the manual package (i.e. manual.HEAD).
Credit
If you are using this code or products in your scientific publication please give a reference to Ofek (2014; ascl.soft 07005).
License
Unless specified otherwise this code and products are released under the GNU general public license version 3.
Installation
See http://weizmann.ac.il/home/eofek/matlab/doc/install.html for installation instruction and additional documentation.
Properties
The SIM object is a strcture-array like object that may contain a large number of images. Each image is stored in a different element of the structure, and each element contains the properties listed below.
The SIM class inherits from the WorldCooSys, HEAD, AstCat, MASK, and ClassPSF classes.
Properties of the image class:
- Im - Image container
- ImageFileName - Image file name container
- BackIm - Background image container
- ErrIm - Error image container
- WeightIm - Weight image container
Properties inherited from the AstCat class:
- Cat - Catalog container
- Col - Structure array of column names and their index
- ColCell - Cell array of column names.
- ColUnits - Cell array of column units.
- SortedBy - Sorted by: column name.
- SortedByCol - Sorted by: Colum index.
- Name
- Source
- Reference
- Version
Properties inherited from the HEAD and WorldCooSys classes:
- Header - A 3 columns cell array of Keyword, Value, Comment
- WCS - Structure of WCS information.
- UserData - User data container
Properties inherited from the MASK class:
- Mask - Mask image container
- MaskDic - Mask dictionary function handle.
Properties inherited from the ClassPSF class:
- PSF - PSF container
- ErrPSF - PSF error container.
- ParPSF - PSF parameters.
- CooPSF - PSF coordinates.
Methods and examples
To see all properties and methods of a SIM object, type "SIM." followed by <tab>.
Constractor
The SIM constractor function can be used to generate an empty SIM object:
% constract a single elemnt SIM object
S=SIM
% constract a 2x1 SIM object
S=SIM(2)
% constrct a 3x2 SIM object
S=SIM(3,2)
Class managment
To check if an object belongs to the SIM class, you can use the isa function or the SIM.issim statcic method:
SIM.issim(1)
S=SIM;
SIM.issim(S)
Note also that the isfield and isstruct functions works on SIM objects.
The following static methods can be used to get the various field names (properties):
SIM.ImageField
SIM.BackField
SIM.ErrField
AstCat.CatField
% etc.
Read images
There are several ways to populate a SIM object with images, or read images into a SIM array.
The first is to populate the fields directly by accessing the fields:
S=SIM;
% populate the SIM object with a scalar image:
S.Im = 1;
% alternatively:
S.(SIM.ImageField) = 1
% populate the SIM image with some images:
S(1).Im = rand(5,5);
S(2).Im = ones(3,3)
The second is to use the getters (under development).
% place holder
The third is to use special functions to pupulate a SIM object. For example, you can use the FITS static class method FITS.read2sim:
S=FITS.read2sim('test_PTF*.fits')
This method have many options to read FITS files (see the FITS class manual for details).
Note that this function, by default, also populates the Header and WCS fields.
Write images
There are several ways to write SIM images into FITS or HDF5 files:
% write a SIM object field into a FITS file
FITS.write(S(1).Im,'FileName.fits');
% delete the file
delete('FileName.fits');
The Header
Since the SIM class inherit from the HEAD class, all the methods available in the HEAD class are available here. Here we provide only a few basic examples (out of the very large number of methods), but a more comprehensive list is available in the manual.HEAD documentation:
% get key
S(1).getkey('DATE-OBS')
% check key values
IsKeyVal=iskeyval(S,'EXPTIME',60,'FILTER','R')
% coordinates conversion:
[RA,Dec]=S(1).xy2coo(rand(5,1),rand(5,1))
% get the JD of the images
julday(S)
Display images
The SIM class interact with the ds9 static class, so SIM images can be easily displayed into the ds9 window. For more details and example see manual.ds9.
Basic operators
Basic binary and unary operators (e.g., +, sin, mean) are overloaded for SIM objects. The basic methods that are being used in the implementation of this overloaded functions are: bfun2sim, ufun2sim, and ufun2scalar. Using these functions it is straightforward to implement and operation on SIM images. Note that these functions allows the user to control on which SIM fields the operators will be executed, the image sub region on which it will be executed, and how to treat the mask images.
Here are some high level examples:
S=FITS.read2sim('test_PTF*.fits')
% add a scalar to each one of the images in the SIM array.
% Note that by default this will be added only to the Im field.
S+1
% while if you like to add 1 to the first SIM element and add 2 to the second:
S+[1;2]
% if you like to do so for the Im and BackIm fields:
Sn=bfun2sim(S,1,@plus,'ExecField',{SIM.ImageField, SIM.BackField});
% to see what happened
S(1).Im(1)
Sn(1).Im(1)
% An example for a unary operator
sin(S(1))
% An example for an operator that returns a scalar
median(S)
% other operators
sin(S).^2 - tan(S) + mean(S)
% You can also execute logical operators on images
S>3
The (almost) full list of supported functions:
plus, minus, power ,power, mtimes, rdivide, ldivide, times, sqrt, uminus, uplus
abs, fix, floor, ceil, transpose, ctranspose, round
case, single, double, int8, int16, int32, uint8, uint16, uint32
gradient, del2
sin, cos, tan, acos, asin, atan2, atan, exp, log10, log
and, eq, ge, gt, le, lt, ne, not, or
angle, circshift, conj, conv2_fft, conv2, conv_fft2, fft2, fft, fftshift, filter, ifft, ifft2, ifftshift, imag, iradon, radon, real,
iqr, prctile, quantile, max, mean, median, min, mode, moment, nanmean, nanmedian, nanstd, rstdnansum, nanvar, skewness, sum, var
Summation of images can be done using coadd method which is further discussed below.
Image scaling can also be done using the scale method (see help for further details).
Bias subtraction
The Bias subtraction process involve identificationof bias images and or overscan regions, removal of bad bias images, constraction of a bias image, identifying bad pixels (based on the bias image), and subtraction of the bias. The SIM/debias method take care of all these steps:
% cd to a directory containing your nightly images
S=FITS.read2sim('*.fits')
S=debias(S,bias(S));
% or
S=debias(S);
% The following example don't trim the images from the overscan region:
S=debias(S2,[],'FinalSec',[]);
Another example, with which fixes a bug in the P200/LFC header keyword:
% debias P200 LFC images
% cd /str1/eran/archive/P200/20111002_LFC/
S=FITS.read2sim('ccd*.0.fits');
S=trim_image(S,[1 1056 1 2063]); % fix bug in Header info
[Sbs,Bias,BiasOS,NIB,sum]=debias(S,[],'IsBiasPar',{'RN',[]},'FinalSec',[1 1024 1 2063]);
See help debias for options and list of output arguments.
Note that the second output argument contains the bias image, its Err image and weight image, as well as the mask image.
The bit mask default dictionary is defined in: MASK.def_bitmask_pipeline. The bias process generate 3 bit masks:
- Bit_BiasNoise0 - Flag pixels with very low noise, typically indicating dead pixels. Default noise is 0.01.
- Bit_BiasNoisy - Flag pixels with very high noise. Default is 4 times the mean of the the std in each pixel.
- Bit_BiasAnom - Anonalous bias level. Default is below 0.7 or above 1.3 of the mean bias level.
Next, we provide some examples related to the lower-level methods that are sometimes useful for more compliacted analysis:
The isbias function can be used to identify bias images and check their validity using several methods. The isbias function first attempt to identify bias images either by header keywords (specifically, the TYPE, IMTYPE, OBSTYPE, IMGTYP, IMGTYPE header keywords), or the file name, or using the image statistics. Then the bias images can be (optionally) validated against a reference bias image and/or the expected readnoise of the detector:
[IsBias,R]=isbias(S);
All the options are documented in the internal help of the function (help SIM/isbias).
The bias method can be used to create the bias image, while the bias_overscan can be used to calculate and remove, line by line, the overscan bias.
Flat correction
There are several functions intended to help with flat field generation and correction. The high level function that do it all is SIM/deflat. This function, can identify flat images, check their validity, create a flat image, identification of bad pixels, and correct all the images by the flat.
Additional low level functions include SIM/isflat to look for flat images and check their validity, and SIM/flat to correct the images by a given flat image.
Example:
% Here Sbs is a SIM array containing the bias subtracted images (see last step)
[Sfc,FlatSim,IsNonFlat,SumFlat]=deflat(Sbs);
The flat correction step also look for bad pixels and update the mask image.
- Bit_FlatNaN - The pixel value is NaN (presuambly due to division by zero).
- Bit_FlatLowNim - The flat field pixel value is based on a small number of images (default is 4).
- Bit_FlatHighStd - The normzlied flat field scatter at this pixel have std larger than some threshold (default is 0.03).
- Bit_FlatLow - The normalized flat field level is lower than some threshold (default is 0.4).
Background and noise image
An important operation is background level and noise estimation. The main function that implement this capability is the background method. This function can use many algorithms. The default algorithm, which is the best based on our experience, is to fit a Gaussian to an histogram of the pixel values.
Example:
S=FITS.read2sim('test_PTF*.fits');
% The following will populate the BackIm and ErrIm felds in the SIM object
S=double(S);
S=background(S);
% an example for some options:
S=background(S,'Block',[256 256]);
% display the back and err images in ds9
ds9(S(1).BackIm,1)
ds9(S(1).ErrIm,2)
Note that by default the background function does not subtract the background (only calcualte it).
You can subtract the background, by either:
Sb=background(S,'SubBack',true);
% or
Sb=sub_background(S);
Source detection and photometry
You can perform source detection and photometry using the mextractor and sextractor methods.
For details see mextractor documentation in manual.mextractor.
A simple example:
S=mextractor(S);
Note that this operation will populate (if does not exist) the BackIm, ErrIm, MASK, PSF and Cat fields. Specifically, the Cat field will be populated with the extracted source catalog.
Cosmic ray identification
There are several approaches for cosmic ray identification - two methods that are supported are edge detection-like algorithms and specifically the van Dokkum (2001) Laplacian Edge Detection, and hypothesis testing. The hypothesis testing approach is fully supported in the mextractor code (see previou section), while the Laplacian Edge Detection is available in the SIM.crdetect_lacos function (see also ImUtil.Im.crdetect_lacos).