TWiki
>
Public Web
>
PacsDocumentation
>
MadMapNotes
(2010-10-04,
BartVandenbussche
)
(raw view)
E
dit
A
ttach
Tags:
create new tag
view all tags
<pre> ====== Notes: ------ HISTORY: 30-Apr-2009: Dave Frayer Initial write-up 14-Oct-2009: Babar Ali Updated API change 12-Jul-2010: Cate Liu Update API documentation 16-Jul-2010: Cate Liu Added new parameter to makeTodArray ====================================== PACS MADmap Users Manual Documentation -------------------------------------- ================= Summary of Method ----------------- MADmap uses a maximum-likelihood technique to build a map from an input Time Order Data (TOD) set by solving a system of linear equations. It is used to remove low-frequency drift ("1/f") noise from bolometer data while preserving the sky signal on large spatial scales. The TOD d = d(t) for t=1...n. is given by (1) d = As+n, where A is the "pointing" matrix which relates each time sample to the corresponding pixel in the sky; A(t,p) for time t=1...n and sky-pixels, p=1...m, s is the sky map s(p), for pixels p=1...m, and n is the noise n(t). The time-time noise covariance matrix N is given by N = <n n^T>, where "^T" represents transpose and "< >" represents the ensemble average. In practice, the noise is assumed to be Gaussian, "piece-wise stationary", and "circulant" which enables the noise properties to be described by its Fourier power spectrum. Under these assumptions, the noise depends only on the time difference between two samples, and N is defined by its first row, N_i0, which is given by the autocorrelation of the noise. (2) N_i0 = <n*n> = <FFT^-1 {[FFT(n)]^2}> where "*" represents convolution and FFT is the Fourier transform. Therefore, the time-time noise covariance matrix N is given by the inverse Fourier transform of its Fourier power spectrum. For Gaussian noise, the probability distribution [P(n)] is given by (ignoring normalization factors): (3) P(n) ~ exp[-0.5 n^T N^-1 n]. Using (1) and (3), the probability of sky (s) given data (d) is: (4) P(s|d) ~ exp[-0.5 (d - As)^T N^-1 (d - As)]. Maximizing this probability yields A^T N^-1 As = A^T N^-1 d, which has a solution: (5) s = (A^T N^-1 A)^-1 A^T N^-1 d. MADmap solves (5) for the sky-map (s) given input TOD (d). MADcap solves the the inverse time-time noise covariance matrix N^-1 (N_i0^-1) without requiring inversion under the assumptions for equation (2) by (6) N_i0^-1 = <FFT^-1 {[FFT(n)]^-2}>. N^-1 effectively acts like a set of convolutions with band-limited kernels. N^-1 is the InvNtt calibration file used by MADmap and is effectively the inverse FFT of 1/(power_spectrum). References: J.-C. Hamilton, 2003, C. R. Physique 4, 871 Ashdown, M. A. J., et al. 2007, A&A, 467, 761 D. Clements, et al. 2006, "SPIRE Mapmaking Algorithm Review Report" ============== Implementation -------------- The implementation is based on three modules: (1) photOffsetCorr (2) photModuleDriftCorrection (3) photGlobalDriftCorrection (4) makeTodArray (5) runMadMap ============== PhotOffsetCorr -------------- I. Summary Removes detector-to-detector offsets before running MADmap. Defaults to a median of input frames cube, but supports removal of an offset calibration file. II. Usage frames = PhotOffsetCorr(Frames frames, int copy, int zeroLevelCorr, Selection regionSelection) frames -- Input frames, MANDATORY, NO default value copy -- default 0: input frames are changed 1: new frames are copied and returned zeroLevelCorr -- default 0: median removal 1: zero level calibration value subtracted regionSelection -- default None: noOp; Otherwise, if zeroLevelCorr=0, do median removal using the given region; if zeroLevelCorr=1, do zero level calibration value subtraction using given region Assumes BAND information is present for selection of proper calibration file. If BAND not present, checks dimensions for Red vs Blue data and defaults to BS for for Blue. III. Functionality frames_OUT = frames_IN - OFFSET where OFFSET is median per detector of frames_IN, or OFFSET is input calibration file. ========================= photModuleDriftCorrection ------------------------- See MADmap User's Guide. ========================= photGlobalDriftCorrection ------------------------- See MADmap User's Guide. ============ makeTodArray ------------ I. Summary Builds time-ordered data (TOD) stream for input into MADmap and derives meta header information of the output skymap. Input data is assumed to be calibrated and flat-fielded. Also prepares the "to's" and "from's" header information for the InvNtt (inverse time-time noise covariance matrix) calibration file. Assumes BAND and OnTarget-flag in the status. II. Usage PacsTodProduct todProd = makeTodArray(Frames frames, Double scale, Double crota2, Boolean optimizedOrientation, Double minRotation, Boolean chunkScanLegs, PacsCal calTree) frames -- Data frames in units of mJy/pixel. Required input meta-data: (1) RA,Dec cubes associated with the frames including the effects of distortion. Assumes this step has been previously done by PhotAssignRaDec, (2) input mask cube which identifies bad pixels, (3) information on band (BS,BL,R), mode (scan/chopped raster), and locations between scan legs for data "chunking". scale -- Default = 1.0. Multiplication pixel scale factor of output sky pixels in relation to nominal PACS detector size (3.2" for blue and 6.4" red). For scale = 1.0, the skymap has square pixels equal to nominal PACS detector size; for scale = 0.5, then 1.6" for blue and 3.2" for red. crota2 -- Default = 0.0 degree. CROTA2 of output skymap. optimizeOrientation -- Default = false when true: do auto map rotation if calculated auto angle > minRotation minRotation -- Default = 15.0 degrees. Minimial angle for auto rotation if optimizeOrientation=true. chunkScanLegs -- Default = true, on-target flags are used to chunk scan legs. i.e. off-target data are not used. calTree -- Default = null, Pacs calibration tree. todProd -- Output product with TOD file name, final output map's Wcs, to/from index for each good data stream segment and the correspondence between tod indices and sky map pixels. The intermediate output TOD file is saved in a directory specified by the property var.hcss.workdir (please note that this property has default value defined in HCSS). It is removed after MADmap finishes and at the exit of HIPE. The body of the TOD file is a byte stream binary data file consisting of header information and TOD data (Reference: http://crd.lbl.gov/~cmc/MADmap/doc/man/MADmap.html). The binary header is four 8byte integers (1) First sample index for TOD data, set to 0. (2) Last sample index for TOD data chunk, set to (n_good_detectors * n_samples) -1. (3) nnObs = Number of detector values per sky pixel during each time sample (for default one-to-one mapping of detectors on to sky pixels, nnObs=1). (4) total number of sky pixels with good data. The binary header is followed by the data in the order of: For each input GOOD detector pixel ("observation"): value (double, 8-byte float) ==v For each sky pixel observed: weight (4-byte float) == w skypixel index (4-byte int) == p (e.g.,.) for good detectors ii=1,nd and time samples kk=1,nt, TOD order is given by: for ii=1,nd for kk=1,nt v[ii,kk] for jj=1,nnObs w[ii,kk] p[ii,kk] Initially for the SPG, we will set nnOBS=1, i.e., use the default one-to-one mapping of input detectors onto sky pixels. The output tod product includes the astrometry of output map using Wcs, in particular meta data keywords such as: CRVAL1 RA Reference position of skymap CRVAL2 Dec Reference position of skymap RADESYS ICRS EQUINOX 2000. CTYPE1 RA---TAN CTYPE2 DEC--TAN CRPIX1 Pixel x value corresponding to CRVAL1 CRPIX2 Pixel y value corresponding to CRVAL2 CDELT1 pixel scale of sky map (=input as default, user parameter) CDELT2 pixel scale of sky map (=input as default, user parameter) CROTA2 PA of image N-axis (=0 as default, user parameter) III. Functionality (1) Build TOD binary data file with format given above. (2) Define the astrometry of output map and save as keywords give above. Default CROTA2=0.0, but if optimizedOrientation=True, then call Maptools to compute the optimial rotation (i.e., for elongated maps). If rotation less than minRotation, then leave map un-rotated with CROTA2=0. (3) BADPIXEL -- For "BAD" dead/bad detector which are always bad, then do not include the detector in TOD calculations Dead/bad detectors will not have good InvNtt data and should be discarded for MADmap. (4) SKYPIX INDICES -- Compute the skypix indices which are derived from the projection of each input pixel onto the output sky grid. The skypix indices are increasing integers representing the location in the sky map with good data. The skypixel indices of the output map must have some data with non-zero weights, must be continuous, must start with 0, and must be sorted with 0 first and the largest index last. (5) GLITCHES -- Set weights=BADwtval for bad data as flagged in the masks (BADPIXEL, SATURATION, GLITCH, UNCLEANCHOP). For (BLIND) Default BADwtval =0.0, but may need to use a small non-zero value (e.g., 0.001) in practice (to avoid MADmap precondition that requires non-zero wts and data for each observed skypixel). Good data should have weights set to 1.0 for initial versions with nnOBS=1. (6) CHUNK per SCAN LEG -- Use the OnTarget status flag to define the boundaries for data chunking per scan leg. The sample start and end boundaries of the TOD indices of each data chunk per detector are needed for the InvNtt headers (the "tos" and "froms"). TOD rules: (1) Chunk for large GAPS (>maxGap) defined by OnTarget (it is assumed that the scan turn-around time will be larger than maxGap, but this is not a requirement). (2) for small GAPS (<=maxGap) defined by OnTarget use the data values for the TOD, but set the TOD weights=BADwtval (i.e., data not effectively used for map, but needed for continuous noise estimate for MADmap). This condition is not expected for properly handled data products upstream, but could exists if there are issues with the pointing products. (3) Include an option to turn-off chunking by scan leg. (7) CHUNK as function of TIME per DETECTOR based on the mask information -- Check TOD stream per detector. For "gaps" of bad data in samples larger than maxGap, then chunked the data. Ignore data streams that have a number of samples less than minGOOD, i.e., each TOD chuck should be larger or equal to minGOOD samples. For gaps, smaller or equal to maxGap, linearly interpolate the TOD values across the gap and set TOD weights=BADwtval (i.e., data not used for map, but needed for continuous noise estimate for MADmap). TOD rules (in this order): (1) throw out initial and end samples that are bad. (2) fill in the small bad GAPS (<=maxGap), weights=BADwtval (3) chuck for large bad GAPS (>maxGap) (4) throw out small chucks (<minGOOD) Initial default maxGap=5 and minGOOD=correlation length of the invntt calibration data. The locations defining the boundaries of the good chunks will need to be stored on a per detector basis (the "tos" and "froms"). Note that (6) and (7) have similar logic [(6) based on OnTarget and (7) based on the mask information]. In both cases, chuck for GAPS>maxGap. For GAPS<=maxGap, linearly interpolate data values across the gaps for (7) [i.e., glitches], but use the data values for (6) [i.e., nothing wrong with the data values]. ========= runMadMap --------- I. Summary Wrapper that runs JAVA MADmap module. MADmap is a maximum likelihood map making routine for converting time ordered data sets into pixelized sky maps. Input TOD data is assumed to be calibrated and flat-fielded and input InvNtt is from calibration tree. (Reference: http://crd.lbl.gov/~cmc/MADmap/doc/man/MADmap.html) II. Usage SimpleImage map = runMadMap(PacsTodProduct tod, PacsCal calTree, Double maxerror, Integer maxiterations, Boolean runNaiveMapper, Boolean useAvgInvntt) Inputs: tod -- PacsTodProduct from makeTodArray, MANDATORY. calTree -- Pacs calibration tree containing calibration InvNtt information stored as an array of size max(n_correlation+1) x n_all_detectors. Each row represents the InvNtt information for each detector. MANDATORY. maxerror -- Default = 1e-5, maximum relative error allowed in PCG routine. maxiterations -- Default = 200, maximum number of iterations in PCG routine. runNaiveMapper -- Default = false, run MadMapper; when true, run NaiveMapper useAvgInvntt -- Default = false, use per-pixel invntt data; when true, use averaged invntt data Output: map -- Output product consisting of following: (1) image -- Sky map image (either a Naive map or MADmap) with Wcs information (2) coverage -- Coverage map corresponding to sky map, with Wcs (3) error -- Uncertainty image associated with the map, with Wcs The error map is currently only made properly for NaiveMapper. Two calls are needed to produce both the MadMap and NaiveMap simple image products (runNaiveMapper=true yields Naivemap product and runNaivMapper=false yields MadMap product). III. Functionality (1) Build InvNtt from input calibration tree. Format of InvNtt chunk: Header-LONG: min_sample starting sample index max_sample last sample index n_correlation correlation width of matrix Data-DOUBLE: invntt(n_correlation+1) The min/max samples are the "tos" and "froms" calculated from a method within makeTodArray. The sample indices need to be consistent between the TOD chunked files and the InvNtt file. (2) Run MADmap module. (3) If runNaiveMapper = true run NaiveMapper module. (4) Put astrometric header information into output products. --- EOF </pre>
E
dit
|
A
ttach
|
Watch
|
P
rint version
|
H
istory
: r1
|
B
acklinks
|
V
iew topic
|
Ra
w
edit
|
M
ore topic actions
Topic revision: r1 - 2010-10-04
-
BartVandenbussche
Public
Log In
Public Web
Create New Topic
Index
Search
Changes
Notifications
Statistics
Preferences
Webs
Public
TWiki