

FUNCTION RESIDUAL_LIMB, img, p, skern, cut_level, prof
;+
; NAME:
;    RESIDUAL_LIMB
;
; PURPOSE:
;    To calculate the symmetric limb darkening function of a fulldisk
;    image given the image.
;
; CALLING SEQUENCE:
;    limb = RESIDUAL_LIMB( img, p, skern, cut_level, prof )
;
; INPUTS:
;    IMG:    Input image.
;    Skern       = maximum size of smoothing kernel
;    Cut_level   = Cutoff level for smoothing filter (%)
;
; OUTPUTS:
;    LIMB:   Limb darkening image.
;    PROF:   1-dimensional limb darkening profile.
;
; KEYWORDS:
;    None.
;
; PROCEDURE:
;    Several tracins are made from the center and out using bilinear
;    interpolation (using 'rad_prof.pro' see BBSO IDL library). These
;    are median filtered for each radius used. The result is smoothed
;    using wider smoothing kernels at the flatter part of the disk.
;    High activity regions (containing fewer than 10 data points in a
;    azimuthal tracin (usually at disc center) are excluded and the
;    center region is instead filled in from a fit to a second order
;    polynomial.
;
; MODIFICATION HISTORY:
;    1997-11-13: Anders Johanneson, Big Bear Solar Observatory
;        Original verion of FD_SYMM.PRO.
;    1999-12-17: Carsten Denker, Big Bear Solar Observatory
;        Intergration in the new BBSO archiving system, variables are passed
;        as structures.
;
;-

;----- radial median filter

r = FIX( ROUND( ( p.xd + p.yd ) / 4. ) )
DIST_CIRCLE, d, [ p.nx, p.ny ], p.xc, p.yc
d = FIX( ROUND( d ) )
index = SORT ( d )
h = HISTOGRAM ( d )
prof = FLTARR( r + 1 )
r1 = 0
FOR i = 0, r DO BEGIN
    r2 = r1 + h( i ) - 1
    prof ( i ) = MEDIAN( img( index( r1 : r2 ) ) )
    r1 = r2 + 1
ENDFOR

;----- fit second order polynomial to inner 25% of the limb darkening profile
;       using 50% instead of 25% seems to give better results

prof = RESIDUAL_DC_INT( prof, r / 4.0, 2 )
;prof = RESIDUAL_DC_INT( prof, r / 2.0, 2 )

;----- filter limb darkening profile

prof = RESIDUAL_SMOOTH( prof, skern, cut_level )

;----- extend profile beyond sun

;========== the next line has modified by Liu Yu, 2002/08/31===huairou==
;limb = REPLICATE( prof( r )/1.15, p.nx, p.ny )
limb = REPLICATE( prof( r )/2., p.nx, p.ny ) ;  for HSOS
r1 = 0
FOR i = 0, r DO BEGIN
    r2 = r1 + h( i ) - 1
    limb( index( r1 : r2 ) ) = prof ( i )
    r1 = r2 + 1
ENDFOR
RETURN, limb

END

