/* @(#)wa_filt.c 17.1.1.1 (ES0-DMD) 01/25/02 17:21:40 */ /*=========================================================================== Copyright (C) 1995 European Southern Observatory (ESO) This program is free software; you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation; either version 2 of the License, or (at your option) any later version. This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details. You should have received a copy of the GNU General Public License along with this program; if not, write to the Free Software Foundation, Inc., 675 Massachusetss Ave, Cambridge, MA 02139, USA. Corresponding concerning ESO-MIDAS should be addressed as follows: Internet e-mail: midas@eso.org Postal address: European Southern Observatory Data Management Division Karl-Schwarzschild-Strasse 2 D 85748 Garching bei Muenchen GERMANY ===========================================================================*/ /****************************************************************************** ** Copyright (C) by European Southern Observatory ******************************************************************************* ** ** UNIT ** ** Version: 17.1.1.1 ** ** Author: Jean-Luc Starck ** ** Date: 02/01/25 ** ** File: wa_filt.c ** ******************************************************************************* ** ** DESCRIPTION filtering by the wavelet transform ** ----------- ** ** ** PARAMETRES ** ---------- ** wavelet transform algoritm (keyword: INPUTI[1]): ** 1 : a trous algorithm with linear scaling function ** The wavelet function is the difference between two resolutions ** 2 : a trous algorithm with B3-spline scaling function ** The wavelet function is the difference between two resolutions ** 3 : algorithm using the FFT without reducing the sampling ** The Fourier transform of the scaling function is a B3-spline ** The wavelet function is the difference between two resolutions ** 4 : pyramidal algorithm with a linear scaling function ** The wavelet function is the difference between two resolutions ** 5 : pyramidal algorithm with a B3-spline scaling function ** The wavelet function is the difference between two resolutions ** 6 : pyramidal algorithm using the FFT: ** The Fourier transform of the scaling function is a B3-spline ** The wavelet function is the difference between two resolutions ** 7 : pyramidal algorithm using the FFT: ** The Fourier transform of the scaling function is a B3-spline ** The wavelet function is the difference between the square of ** two resolutions ** 8 : Mallat's Algorithm with biorthogonal filters ** ** Type_Filter: (keyword: INPUTI[2]): defines the choosen filtering ** FILTER_TRESHOLD 1 ==> Thresholding ** FILTER_HIERARCHICAL_TRESHOLD 2 ==> Hierarchical Thresholding ** FILTER_HIERARCHICAL 3 ==> Hierarchical Wiener Filtering ** FILTER_MULTI_RES_WIENER 4 ==> Multiresolution Wiener Filtering ** ** The best results are obtained with the Hierarchical Wiener Filtering ** ** N_Sigma: (keyword: INPUTR[1]): If we threshold, the level is estimated by: ** Level = N_Sigma * Noise_Standart_Deviation ** N_Sigma = 3 is a standart value ** ** Nbr_Iter: (keyword: INPUTI[3]): ** if we threshold, the reconstruction is iterative ** Nbr_Iter between 6 and 10 is generally good ** Nbr_Iter = 1 ==> no iterative reconstruction ** ** Nbr_Plan: (keyword: INPUTI[4]): number of scales ** it is not necessary to take values for Nbr_Plan superior to 5 ** ** Sigma_Noise: (keyword: INPUTR[3]): ** standard deviation of the noise ** if Sigma_Noise = 0 then the standard deviation of the noise ** is estimated automatically ** ** RESULTS a file is created which contains the filtered image ** ------- ** ******************************************************************************/ static char sccsid[] = "@(#)wa_filt.c 17.1.1.1 02/01/25 ESO @(#)"; #include #include #include #include "Def_Math.h" #include "Def_Mem.h" #include "Def_Wavelet.h" #define VISU_DATA FALSE /***************************************************************************/ main() { string File_Name_Imag, File_Name_Imag_Out; int Type_Transform, Nbr_Plan; float *Imag; int Nl, Nc, Type_Filter, Nbr_Iter = 1; float Fc = 0.5, N_Sigma, Noise; int Unit; int Stat, Null, Actvals, Buffer_Int[4], Maxvals, Felem; double Buffer_Float; /* Initialisation */ SCSPRO("transform"); /* read the image file name */ Felem = 1; Maxvals = 60; Stat = SCKGETC("IN_A", Felem, Maxvals, &Actvals, File_Name_Imag); /* read the output file name */ Stat = SCKGETC("OUT_A", Felem, Maxvals, &Actvals, File_Name_Imag_Out); /* read the transform type, chosen filtering, number of iterations and the number of scales */ Felem = 1; Maxvals = 4; Stat = SCKRDI("INPUTI",Felem, Maxvals, &Actvals, Buffer_Int, &Unit, &Null); Type_Transform = Buffer_Int[0]; Type_Filter = Buffer_Int[1]; Nbr_Iter = Buffer_Int[2]; Nbr_Plan = Buffer_Int[3]; /* read the N_Sigma parameter and noise standard deviation */ Felem = 1; Maxvals = 1; Stat = SCKRDR("INPUTR", Felem, Maxvals, &Actvals, &N_Sigma, &Unit, &Null); Felem = 2; Stat = SCKRDR("INPUTR", Felem, Maxvals, &Actvals, &Noise, &Unit, &Null); /* load the input image */ io_read_file_to_pict_f (File_Name_Imag, &Imag, &Nl, &Nc); #if VISU_DATA printf ("Type_Filter = %d\n", Type_Filter); printf ("N_Sigma = %f\n", N_Sigma); printf ("Type_Transform = %d\n", Type_Transform); printf ("Nbr_Plan = %d\n", Nbr_Plan); #endif /* filtering*/ wave_filter_imag (Imag, Nl, Nc, Imag, N_Sigma, Type_Filter, Nbr_Iter, Type_Transform, Nbr_Plan, Fc, Noise); /* save the result */ io_write_pict_f_to_file (File_Name_Imag_Out, Imag, Nl, Nc); free ((char *) Imag); /* End */ SCSEPI(); } /***************************************************************************/