mxlib
c++ tools for analyzing astronomical data and other tasks by Jared R. Males. [git repo]
|
An implementation of the Karhunen-Loeve Image Processing (KLIP) algorithm.
KLIP[24] is a principle components analysis (PCA) based technique for PSF estimation.
_realT | is the floating point type in which to do calculations |
_derotFunctObj | the ADIobservation derotator class. |
_evCalcT | the real type in which to do eigen-decomposition. Should generally be double for stable results. |
Definition at line 94 of file KLIPreduction.hpp.
#include <improc/KLIPreduction.hpp>
Public Member Functions | |
KLIPreduction () | |
Default c'tor. More... | |
KLIPreduction (const std::string &dir, const std::string &prefix, const std::string &ext) | |
Construct and load the target file list. More... | |
KLIPreduction (const std::string &fileListFile) | |
Construct using a file containing the target file list. More... | |
KLIPreduction (const std::string &dir, const std::string &prefix, const std::string &ext, const std::string &RDIdir, const std::string &RDIprefix, const std::string &RDIext="") | |
KLIPreduction (const std::string &fileListFile, const std::string &RDIfileListFile) | |
Construct using a file containing the target file list and a file containing the RDI target file list. More... | |
void | meanSubtract (eigenCube< realT > &rims, eigenCube< realT > &tims, std::vector< realT > &sds) |
Subtract the basis mean from each of the images. More... | |
int | regions (std::vector< realT > minr, std::vector< realT > maxr, std::vector< realT > minq, std::vector< realT > maxq) |
Run KLIP in a set of geometric search regions. More... | |
int | regions (realT minr, realT maxr, realT minq, realT maxq) |
Run KLIP in a geometric search region. More... | |
void | worker (eigenCube< realT > &rims, eigenCube< realT > &tims, std::vector< size_t > &idx, realT dang, realT dangMax) |
Public Member Functions inherited from mx::improc::ADIobservation< _realT, _derotFunctObj > | |
ADIobservation (const std::string &dir, const std::string &prefix, const std::string &ext) | |
ADIobservation (const std::string &fileListFile) | |
ADIobservation (const std::string &dir, const std::string &prefix, const std::string &ext, const std::string &RDIdir, const std::string &RDIprefix, const std::string &RDIext="") | |
ADIobservation (const std::string &fileListFile, const std::string &RDIfileListFile) | |
int | readFiles () |
Read in the target files. More... | |
virtual int | postReadFiles () |
Post target read actions, including fake injection. More... | |
virtual int | postCoadd () |
Post target coadd actions. More... | |
int | readRDIFiles () |
Read in the RDI files. More... | |
virtual int | postRDIReadFiles () |
Post reference read actions, including fake injection. More... | |
virtual int | postRDICoadd () |
Post reference coadd actions. More... | |
int | readPSFSub (const std::string &dir, const std::string &prefix, const std::string &ext, size_t nReductions) |
Read in already PSF-subtracted files. More... | |
virtual void | makeMaskCube () |
Populate the mask cube which is used for post-processing. Derived classes can do this as appropriate, e.g. by rotating the mask. More... | |
void | derotate () |
De-rotate the PSF subtracted images. More... | |
int | injectFake (eigenCube< realT > &ims, std::vector< std::string > &fileList, derotFunctObj &derotF, realT RDIfluxScale, realT RDISepScale) |
Inect the fake plants. More... | |
int | injectFake (eigenImageT &fakePSF, eigenCube< realT > &ims, int image_i, realT derotAngle, realT PA, realT sep, realT contrast, realT scale, realT RDIfluxScale, realT RDISepScale) |
Public Member Functions inherited from mx::improc::HCIobservation< _realT > | |
int | readPSFSub (const std::string &dir, const std::string &prefix, const std::string &ext, size_t nReductions) |
Read in already PSF-subtracted files. More... | |
HCIobservation () | |
Default c'tor. More... | |
HCIobservation (const std::string &dir, const std::string &prefix, const std::string &ext) | |
Construct and load the target file list. More... | |
HCIobservation (const std::string &fileListFile) | |
Construct using a file containing the target file list. More... | |
HCIobservation (const std::string &dir, const std::string &prefix, const std::string &ext, const std::string &RDIdir, const std::string &RDIprefix, const std::string &RDIext="") | |
Construct and load the target file list and the RDI file list. More... | |
HCIobservation (const std::string &fileListFile, const std::string &RDIfileListFile) | |
Construct using a file containing the target file list and a file containing the RDI target file list. More... | |
void | load_fileList (const std::string &dir, const std::string &prefix, const std::string &ext) |
Load the file list. More... | |
void | load_fileList (const std::string &fileListFile) |
Load the file list from a file. More... | |
void | load_RDIfileList (const std::string &dir, const std::string &prefix, const std::string &ext) |
Load the RDI basis file list. More... | |
void | load_RDIfileList (const std::string &fileListFile) |
Load the RDI basis file list from a file. More... | |
int | threshold (std::vector< std::string > &fileList, const std::string &qualityFile, realT qualityThreshold) |
Read the image qualities from a qualityFile and apply the threshold to a fileList. More... | |
int | readFiles () |
Read the list of files, cut to size, and preprocess. More... | |
int | readRDIFiles () |
Read the list of reference files, cut to size, and preprocess. More... | |
void | coaddImages (int coaddCombineMethod, int coaddMaxImno, int coaddMaxTime, std::vector< std::string > &coaddKeywords, std::vector< double > &imageMJD, std::vector< fits::fitsHeader > &heads, eigenCube< realT > &ims) |
Coadd the images. More... | |
void | readMask () |
Read the mask file, resizing to imSize if needed. More... | |
void | preProcess (eigenCube< realT > &ims) |
Do the pre-processing. More... | |
int | readWeights () |
Read the image weights from m_weightFile. More... | |
void | combineFinim () |
Combine the images into a single final image. More... | |
void | outputPreProcessed () |
Output the pre-processed target images. More... | |
void | outputRDIPreProcessed () |
Output the pre-processed reference images. More... | |
void | stdFitsHeader (fits::fitsHeader &head) |
Fill in the HCIobservation standard FITS header. More... | |
void | writeFinim (fits::fitsHeader *addHead=0) |
Write the final combined image to disk. More... | |
void | outputPSFSub (fits::fitsHeader *addHead=0) |
Write the PSF subtracted images to disk. More... | |
Public Attributes | |
int | m_meanSubMethod {HCI::imageMean} |
Specify how the data are centered for PCA within each search region. More... | |
std::vector< int > | m_Nmodes |
Specifies the number of modes to include in the PSF. More... | |
realT | m_minDPx {0} |
Specify the minimum pixel difference at the inner edge of the search region. More... | |
realT | m_maxDPx {0} |
Specify the maximum pixel difference at the inner edge of the search region. More... | |
int | m_excludeMethod {HCI::excludeNone} |
Controls how reference images are excluded, if at all, from the covariance matrix for each target image based on a minimum criterion. More... | |
int | m_excludeMethodMax {HCI::excludeNone} |
Controls how reference images are excluded, if at all, from the covariance matrix for each target image based on a maximum criterion. More... | |
int | m_includeRefNum {0} |
Number of reference images to include in the covariance matrix. More... | |
int | m_includeMethod {HCI::includeAll} |
Controls how number of included images is calculated. More... | |
Public Attributes inherited from mx::improc::ADIobservation< _realT, _derotFunctObj > | |
int | m_fakeMethod {HCI::single} |
Method for reading fake files, either HCI::single or HCI::list. More... | |
std::string | m_fakeFileName |
FITS file containing the fake planet PSF to inject or a list of fake images. More... | |
std::string | m_fakeScaleFileName |
One-column text file containing a scale factor for each point in time. More... | |
std::vector< realT > | m_fakeSep |
Separation(s) of the fake planet(s) More... | |
std::vector< realT > | m_fakePA |
Position angles(s) of the fake planet(s) More... | |
std::vector< realT > | m_fakeContrast |
Contrast(s) of the fake planet(s) More... | |
realT | m_RDIFluxScale {1} |
Flux scaling to apply to fake planets injected in RDI. Would depend on the assumed spectrum in SDI. More... | |
realT | m_RDISepScale {1} |
Scaling to apply to fake planet separation in RDI. Would be ratio of wavelengths for SDI. More... | |
Public Attributes inherited from mx::improc::HCIobservation< _realT > | |
eigenCube< realT > | m_tgtIms |
The target image cube. More... | |
int | m_Nims {0} |
Number of images in m_tgtIms. More... | |
int | m_Nrows {0} |
Number of rows of the images in m_tgtIms. More... | |
int | m_Ncols {0} |
Number of columns of the images in m_tgtIms. More... | |
int | m_Npix {0} |
Pixels per image, that is Nrows*Ncols. More... | |
std::vector< double > | m_imageMJD |
Vector of target image times, in MJD. More... | |
std::vector< fits::fitsHeader > | m_heads |
Vector of FITS headers, one per file, populated with the values for the keywords. More... | |
bool | m_filesRead {false} |
Whether or not the m_fileList has been read. More... | |
bool | m_filesDeleted {false} |
Whether or not the specified files have been deleted from m_fileList. More... | |
std::vector< realT > | m_comboWeights |
Vector to hold the image weights read from the m_weightFile. More... | |
eigenCube< realT > | m_refIms |
The optional reference image cube. More... | |
std::vector< double > | m_RDIimageMJD |
Vector of reference image times, in MJD. More... | |
std::vector< fits::fitsHeader > | m_RDIheads |
Vector of FITS headers, one per reference file, populated with the values for the keywords. More... | |
bool | m_RDIfilesRead {false} |
Whether or not the reference files have been read. More... | |
bool | m_RDIfilesDeleted {false} |
Whether or not the specified files have been deleted from m_RDIfileList. More... | |
std::vector< eigenCube< realT > > | m_psfsub |
The PSF subtracted images. More... | |
eigenCube< realT > | m_finim |
The final combined images, one for each cube in psfsub. More... | |
std::vector< std::string > | m_fileList |
The list of files to read in. More... | |
int | m_deleteFront {0} |
Specify how many files from m_fileList to delete from the front of the list. More... | |
int | m_deleteBack {0} |
Specify how many files from m_fileList to delete from the back of the list. More... | |
std::string | m_qualityFile |
File containing 2 space-delimited columns of fileVame qualityValue pairs. More... | |
realT | m_qualityThreshold {0} |
Threshold to apply to qualityValues read from qualityFile. More... | |
bool | m_thresholdOnly {false} |
Just prints the names and qualities of the files which pass threshold, and stop. More... | |
std::string | m_MJDKeyword {"DATE-OBS"} |
Name of the keyword to use for the image date. More... | |
bool | m_MJDisISO8601 {true} |
Whether or not the date is in ISO 8601 format. More... | |
realT | m_MJDUnits {1.0} |
If the date is not ISO 8601, this specifies the conversion to Julian Days (i.e. seconds to days) More... | |
std::vector< std::string > | m_keywords |
Vector of FITS header keywords to read from the files in m_fileList. More... | |
int | m_imSize {0} |
Set the image size. Images are cut down to this size after reading. More... | |
std::vector< std::string > | m_RDIfileList |
The list of files to read in for the RDI basis. More... | |
int | m_RDIdeleteFront {0} |
Specify how many files from m_RDIfileList to delete from the front of the list. More... | |
int | m_RDIdeleteBack {0} |
Specify how many files from m_RDIfileList to delete from the back of the list. More... | |
std::string | m_RDIqualityFile |
File containing 2 space-delimited columns of fileMame qualityValue pairs for the reference images. More... | |
realT | m_RDIqualityThreshold {0} |
Threshold to apply to qualityValues read from qualityFile. More... | |
std::vector< std::string > | m_RDIkeywords |
Vector of FITS header keywords to read from the files in m_fileList. More... | |
int | m_coaddCombineMethod {HCI::noCombine} |
Determine how to coadd the raw images. More... | |
int | m_coaddMaxImno {0} |
Maximum number of images to coadd. More... | |
realT | m_coaddMaxTime {0} |
Maximum elapsed time over which to coadd the images. More... | |
std::vector< std::string > | m_coaddKeywords |
The values of these keywords will be averaged and replaced. More... | |
std::string | m_maskFile |
Specify a mask file to apply. More... | |
eigenImageT | m_mask |
The mask. More... | |
eigenCube< realT > | m_maskCube |
A cube of masks, one for each input image, which may be modified versions (e.g. rotated) of mask. More... | |
bool | m_skipPreProcess {false} |
Don't do any of the pre-processing steps (including coadding). More... | |
bool | m_preProcess_beforeCoadd {false} |
controls whether pre-processing takes place before or after coadding More... | |
bool | m_preProcess_mask {true} |
If true, the mask is applied during each pre-processing step. More... | |
bool | m_preProcess_subradprof {false} |
If true, a radial profile is subtracted from each image. More... | |
realT | m_preProcess_azUSM_azW {0} |
Azimuthal boxcar width for azimuthal unsharp mask [pixels]. More... | |
realT | m_preProcess_azUSM_maxAz {45} |
Mazimum azimuthal boxcar width for azimuthal unsharp mask [degrees]. More... | |
realT | m_preProcess_azUSM_radW {0} |
Radial boxcar width for azimuthal unsharp mask [pixels]. More... | |
int | m_preProcess_medianUSM_fwhm {0} |
Kernel FWHM for symmetric box median unsharp mask (USM) More... | |
realT | m_preProcess_gaussUSM_fwhm {0} |
Kernel FWHM for symmetric Gaussian unsharp mask (USM) More... | |
std::string | m_preProcess_outputPrefix |
Set path and file prefix to output the pre-processed images. More... | |
bool | m_preProcess_only {false} |
If true, then we stop after pre-processing. More... | |
int | m_combineMethod {HCI::meanCombine} |
Determine how to combine the PSF subtracted images. More... | |
std::string | m_weightFile |
Specifies a file containing the image weights, for combining with weighted mean. More... | |
realT | m_sigmaThreshold {0} |
The standard deviation threshold used if combineMethod == HCI::sigmaMeanCombine. More... | |
realT | m_minGoodFract {0.0} |
The minimum fraction of good (un-masked) pixels to include in the final combination (0.0 to 1.0). If not met, then the pixel will be NaN-ed. More... | |
int | m_doWriteFinim {1} |
Set whether the final combined image is written to disk. More... | |
std::string | m_outputDir |
The directory where to write output files. More... | |
std::string | m_finimName {"finim_"} |
The base file name of the output final image. More... | |
bool | m_exactFinimName {false} |
Use m_finimName exactly as specified, without appending a number or an extension. More... | |
bool | m_doOutputPSFSub {false} |
Controls whether or not the individual PSF subtracted images are written to disk. More... | |
std::string | m_PSFSubPrefix |
Prefix of the FITS file names used to write individual PSF subtracted images to disk if m_doOutputPSFSub is true. More... | |
Additional Inherited Members | |
Public Types inherited from mx::improc::HCIobservation< _realT > | |
typedef _realT | realT |
The arithmetic type used for calculations. Does not have to match the type in images on disk. More... | |
typedef Eigen::Array< realT, Eigen::Dynamic, Eigen::Dynamic > | eigenImageT |
The Eigen image array type basted on realT. More... | |
mx::improc::KLIPreduction< _realT, _derotFunctObj, _evCalcT >::KLIPreduction |
Default c'tor.
Definition at line 321 of file KLIPreduction.hpp.
mx::improc::KLIPreduction< _realT, _derotFunctObj, _evCalcT >::KLIPreduction | ( | const std::string & | dir, |
const std::string & | prefix, | ||
const std::string & | ext | ||
) |
Construct and load the target file list.
Populates the m_fileList vector by searching on disk for files which match "dir/prefix*.ext". See load_fileList
[in] | dir | the directory to search for files. |
[in] | prefix | the prefix of the files |
[in] | ext | the extension of the files |
Definition at line 326 of file KLIPreduction.hpp.
|
explicit |
Construct using a file containing the target file list.
Populates the m_fileList vector by reading the file, which should be a single column of new-line delimited file names.
[in] | fileListFile | the full path to a file containing a list of files |
Definition at line 334 of file KLIPreduction.hpp.
mx::improc::KLIPreduction< _realT, _derotFunctObj, _evCalcT >::KLIPreduction | ( | const std::string & | dir, |
const std::string & | prefix, | ||
const std::string & | ext, | ||
const std::string & | RDIdir, | ||
const std::string & | RDIprefix, | ||
const std::string & | RDIext = "" |
||
) |
Populates the m_fileList vector by searching on disk for files which match "dir/prefix*.ext". See load_fileList
Populates the m_RDIfileList vector by searching on disk for files which match "RDIdir/RDIprefix*.RDIext". See load_RDIfileList
[in] | dir | the directory to search for files. |
[in] | prefix | the prefix of the files |
[in] | ext | the extension of the files, default is .fits |
[in] | RDIdir | the directory to search for the reference files. |
[in] | RDIprefix | the initial part of the file name for the reference files. Can be empty "". |
[in] | RDIext | [optional] the extension to append to the RDI file name, must include the '.'. If empty "" then same extension as target files is used. |
Definition at line 340 of file KLIPreduction.hpp.
|
explicit |
Construct using a file containing the target file list and a file containing the RDI target file list.
Populates the m_fileList vector by reading the file, which should be a single column of new-line delimited file names.
Populates the m_RDIfileList vector by reading the file, which should be a single column of new-line delimited file names.
[in] | fileListFile | a file name path to read for the target file names. |
[in] | RDIfileListFile | a file name path to read for the reference file names. |
Definition at line 352 of file KLIPreduction.hpp.
void mx::improc::KLIPreduction< _realT, _derotFunctObj, _evCalcT >::meanSubtract | ( | eigenCube< realT > & | rims, |
eigenCube< realT > & | tims, | ||
std::vector< realT > & | sds | ||
) |
Subtract the basis mean from each of the images.
The mean is subtracted according to m_meanSubMethod.
rims | [in.out] The reference images. These are mean subtracted on output. | |
tims | [in.out] The target images, which can be the same cube as rims (tested by pointer comparison), in which case they will be ignored. Mean subtractedon output. | |
[out] | sds | The standard deviation of the mean subtracted refernce images. |
Definition at line 364 of file KLIPreduction.hpp.
References mx::improc::eigenCube< dataT >::image(), mx::improc::HCI::imageMean, mx::improc::HCI::imageMedian, mx::improc::eigenCube< dataT >::mean(), mx::improc::HCI::meanImage, mx::improc::eigenCube< dataT >::median(), and mx::improc::HCI::medianImage.
|
inline |
Run KLIP in a geometric search region.
This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts.
[in] | minr | |
[in] | maxr | |
[in] | minq | |
[in] | maxq |
Definition at line 262 of file KLIPreduction.hpp.
References mx::improc::KLIPreduction< _realT, _derotFunctObj, _evCalcT >::regions().
int mx::improc::KLIPreduction< _realT, _derotFunctObj, _evCalcT >::regions | ( | std::vector< realT > | minr, |
std::vector< realT > | maxr, | ||
std::vector< realT > | minq, | ||
std::vector< realT > | maxq | ||
) |
Run KLIP in a set of geometric search regions.
The arguments are 4 vectors, where each entry defines one component of the search region.
[in] | minr | |
[in] | maxr | |
[in] | minq | |
[in] | maxq |
Definition at line 447 of file KLIPreduction.hpp.
References mx::improc::cutImageRegion(), mx::math::dtor(), mx::improc::HCI::excludeAngle, mx::improc::HCI::excludeImno, mx::improc::HCI::excludeNone, mx::improc::HCI::excludePixel, mx::sys::get_curr_time(), mx::improc::eigenCube< dataT >::image(), and mx::fits::fitsFile< dataT >::write().
Referenced by mx::improc::KLIPreduction< _realT, _derotFunctObj, _evCalcT >::regions().
void mx::improc::KLIPreduction< _realT, _derotFunctObj, _evCalcT >::worker | ( | eigenCube< realT > & | rims, |
eigenCube< realT > & | tims, | ||
std::vector< size_t > & | idx, | ||
realT | dang, | ||
realT | dangMax | ||
) |
omp_get_num_threads();
Definition at line 780 of file KLIPreduction.hpp.
References mx::improc::eigenCube< dataT >::asVectors(), mx::math::calcKLModes(), mx::improc::eigenCube< dataT >::cube(), mx::math::eigenSYRK(), mx::improc::HCI::excludeNone, mx::sys::get_curr_time(), mx::ipc::ompLoopWatcher< _outputT, _printPretty, _printLoops, _printPercent, _printNLine, _time >::incrementAndOutputStatus(), mx::improc::insertImageRegion(), and mx::fits::fitsFile< dataT >::write().
int mx::improc::KLIPreduction< _realT, _derotFunctObj, _evCalcT >::m_excludeMethod {HCI::excludeNone} |
Controls how reference images are excluded, if at all, from the covariance matrix for each target image based on a minimum criterion.
Can have the following values:
Definition at line 197 of file KLIPreduction.hpp.
int mx::improc::KLIPreduction< _realT, _derotFunctObj, _evCalcT >::m_excludeMethodMax {HCI::excludeNone} |
Controls how reference images are excluded, if at all, from the covariance matrix for each target image based on a maximum criterion.
Can have the following values:
Definition at line 206 of file KLIPreduction.hpp.
int mx::improc::KLIPreduction< _realT, _derotFunctObj, _evCalcT >::m_includeMethod {HCI::includeAll} |
Controls how number of included images is calculated.
The number of included images is calculated after exclusion is complete. Can have the following values:
Definition at line 224 of file KLIPreduction.hpp.
int mx::improc::KLIPreduction< _realT, _derotFunctObj, _evCalcT >::m_includeRefNum {0} |
Number of reference images to include in the covariance matrix.
If > 0, then at most this many images, determined by highest cross-correlation, are included. This is determined after rotational/image-number exclusion. If == 0, then all reference images are included.
Definition at line 213 of file KLIPreduction.hpp.
realT mx::improc::KLIPreduction< _realT, _derotFunctObj, _evCalcT >::m_maxDPx {0} |
Specify the maximum pixel difference at the inner edge of the search region.
Definition at line 188 of file KLIPreduction.hpp.
int mx::improc::KLIPreduction< _realT, _derotFunctObj, _evCalcT >::m_meanSubMethod {HCI::imageMean} |
Specify how the data are centered for PCA within each search region.
Can have the following values:
Definition at line 164 of file KLIPreduction.hpp.
realT mx::improc::KLIPreduction< _realT, _derotFunctObj, _evCalcT >::m_minDPx {0} |
Specify the minimum pixel difference at the inner edge of the search region.
Definition at line 185 of file KLIPreduction.hpp.
std::vector<int> mx::improc::KLIPreduction< _realT, _derotFunctObj, _evCalcT >::m_Nmodes |
Specifies the number of modes to include in the PSF.
The output image is a cube with a plane for each entry in m_Nmodes. Only the number of eigenvalues required for the maximum value of m_Nmodes are calculated, so this can have an important impact on speed.
Can be initialized as:
Definition at line 180 of file KLIPreduction.hpp.