$SSW_HESSI/idl/image/hsi_mem_sato.pro
IDL> o=hsi_image() IDL> o->set,image_alg='mem_sato' IDL> image_mem_sato=o->getdata()
image_mem_sato will contain a two-dimensional array of the reconstructed cartesian (rectangular) image with image dimensions, pixel sizes, and location specified inside the image object.
Sato, Kosugi, and Makishima 1999 PASJ, 51, 127, Improvement of YOHKOH Hard X-Ray Imaging
Maximum Entropy Method - Jim McTiernan (16-Oct-1999)
Comparison of HXT Fortran & IDL routines - Jim McTiernan (25 June 1998)
MEM_sato is a maximum entropy reconstruction of an image that works directly with counts (from the calibrated event list). The algorithm fits the distribution of brightness in the image and the total image brightness (btot) separately. It is based on the Yohkoh Hard X-ray Telescope (HXT) reconstruction algorithm developed by Jun Sato.
If the time interval and energy range are such that the counts are only a few per time bin for any detector then chi squared for the best fit may not be 1. Over-fitting of the noise on the counts can result from this and so give a broken-up image.
Currently, χ2 is used instead of the Cash Statistic when the number of counts in a time bin is very low and the Poisson probability distribution no longer approximates to a Gaussian.
All mem_sato control parameters are prefixed by sato_ and have a similar meaning to those of mem_vis, prefixed by vis_. However, the absolute value of the parameters may have a different effect. A full listing is given in the GUI Guide. The principal parameters are as follows:
Parameter | Function | Default |
sato_lnorm | Controls the initial strength of the entropy constraint - the value of lambda. A small value corresponds to a strong smoothness constraint, and a large value corresponds to a weak constraint. A large value allows for fast reconstruction of an image containing point sources, but will result in very poor broken up images if the sources are extended. | 0.1 |
sato_lambda_max | During image reconstruction, each lambda iteration corresponds to a weakening of the smoothness constraint, allowing a better fit (smaller χ2) to the data. This parameter sets the maximum value that the lambda iteration can reach at which point mem_sato will exit and return the last image. Note that the lambda iteration step is not fixed, it accelerates as the image reconstruction progresses. | 20 |
sato_iter_max | Maximum number of iterations before procedure stops. | 30 |
sato_chi_limit | Sets the χ2 at which the image reconstruction will cease. | 1.1 |
sato_no_chi2 | If set to 1, mem_sato iterates to sato_lambda_max, and ignores sato_chi_limit. | 0 |
Note: mem_sato may terminate before either the chi_limit or lambda_max is reached if χ2 increases or persistently decreases too slowly.
None
mem_sato opens a plot window containing the reconstructed image at the end of the most recent lambda iteration. The image is cartesian (rectangular) with its center at the xy-offset from Sun-center specified in the image object. Pixel sizes and image dimensions are also set according to the values in the object.
The IDL command line or console window will display the current lambda iteration number and the iteration within that lambda value together with its χ2 value. It also shows the current value of btot - the total image brightness.
The reconstructed image is returned by o->getdata() (variable d in the example).
Beware of images that have not reached χ2 of ~1 or that have a "broken up" or "patterned" appearance. The following remedies are recommended:
If the final image has the reduced χ2 > 1, try increasing sato_lambda_max if it was reached, and running the reconstruction again. A larger value of lnorm may speed up the reconstruction, but at the price of a broken up image. If the iterations stopped because χ2 wasn't decreasing fast enough, or was increasing, try altering the value of sato_lnorm.
If the final image appears broken up or shows patterns, try reducing sato_lnorm by a factor of ten and running the reconstruction again.
Use the Cash Statistic rather than χ2 since it is more appropriate when the number of counts in a time bin is very low and the Poisson probability distribution no longer approximates to a Gaussian.