l1a.c

Go to the documentation of this file.

/*
 
$Header: /app/shared/RCS/irix-5.2/seawifsd/src/hdfio/L1A.V4.5/l1a.c,v 4.18 1996/01/17 16:53:52 seawifsd Exp seawifsd $
$Log: l1a.c,v $
Revision 4.18  1996/01/17 16:53:52  seawifsd
added ON_TIME()/OFF_TIME() to report execution elapse time.
 
Revision 4.17  1995/12/19 14:51:04  seawifsd
fixed bug in testing main program. The bug only occurs when
L1B HDF file was not created(specified). Fixes are done to
use the fid of L1A HDF file instead of that of L1B.
 
Revision 4.16  1995/12/07 19:21:01  seawifsd
dded support to compile with calibration and straylight correction routines.
added snode,enode,startclat,startclon,endclat,endclon and
made it IO_SPEC_V4.4 and PROD_SPEC_V2.8 compliant.
added l2_flags, checked and set HIGHLT1 flag. Created additional
SDS to store L1B data in the test main program.
 
Revision 4.15  1995/05/12 13:53:00  seawifsd
1. code conformed to specifications IO_SPEC_V43, PROD_SPEC_V27,
and NONPROD_SPEC_V12(no funtionality changes, just eliminate
code related to versions before IO_SPEC_V42)
 
Revision 4.14  1995/05/04 15:28:24  seawifsd
changed all references to 'exist.h' to 'generic.h'.
 
Revision 4.13  1995/04/03 16:06:36  seawifsd
modified testing part of the l1a.c for more consistent and accurate output.
 
Revision 4.12  1995/02/15 18:51:23  seawifsd
added additional code to implement the objects 'scsol_z', 'entry_year',
'entry_day', and 'csol_z'.
renamed macro  LAKSHMI to EXTRA_FOR_BROWSE.
specifically declared the datatype of fid[] and fid2[].
added an optional third command parameter <skip> in the testing
main program to indicate number of scan line to skip in displaying
l1a science data. Corresponding print statements were revised to
add this support and the support of TRANSPOSE_IMAGE macro.
 
Revision 4.11  1995/01/27 18:20:33  seawifsd
changed the order of mirror from mirror[bands][sides] to mirror[sides][bands]
 
Revision 4.10  1995/01/17 20:03:18  seawifsd
Jan. 17, 1994 V4.1
 
Revision 4.2  1995/01/17 14:34:50  seawifsd
added calling FORTRAN cdata_() before geonav_() because the changes
in MOPS geonav.f code.
 
Revision 4.1  1995/01/17 14:15:37  seawifsd
Jan. 9, 1994, 4.0
 
Revision 3.6  1994/12/15 16:06:42  seawifsd
code conformed to IO_SPEC_V41, PROD_SPEC_V24 and NONPROD_SPEC_V10(no
changes were needed from IO_SPEC_V40)
 
Revision 3.5  1994/12/15 16:04:57  seawifsd
made sure that the pointers passed into the i/o routines in the test
main program are all null by declaring as static.
 
Revision 3.4  1994/12/06 19:31:41  seawifsd
made the code comformed to IO_SPEC_V40 and PROD_SPEC_V24.
elliminated all references to IO_SPEC_V10 and related features under that
version.
made the L1A image data to be pixel interlaced(to fit the spec.) instead
of the original band interlaced.(by define TRANSPOSE_IMAGE in wifs_conf.h)
 
Revision 3.3  1994/11/08 18:45:55  seawifsd
Nov. 8, 1994, 3.3a3
 
Revision 3.3  1994/11/08 15:03:53  seawifsd
Nov. 8, 1994, 3.3a2
 
Revision 1.1.1.5  1994/11/03 20:24:00  frank
renamed MIN macro to MINV.
 
Revision 1.1.1.4  1994/11/01 16:02:28  frank
added processing of 'meta_l1a.infiles' because of the new requirement
from L1A Browse routines.
 
Revision 1.1.1.3  1994/10/21 19:43:51  frank
turned off some message printing by using PRINTF instead of printf.
 
Revision 1.1.1.2  1994/10/05 20:14:01  frank
1. updated L1A source code to fit the specification of the
"SeaWiFS OPERATIONAL ARCHIVE PRODUCT SPECIFICATIONS" v 1.0, 07/22/94,
and the "INTERFACE SPECIFICATIONS FOR SeaWiFS OPERATIONAL PRODUCT
INPUT, OUTPUT, AND SENSOR CALIBRATION SOFTWARE" v 3.3 07/12/94
 
Revision 1.1.1.1  1994/05/23 19:04:43  frank
variables ptr_orb_vec, ptr_l_vert, ptr_sun_ref, ptr_att_ang, ptr_sen_mat,
ptr_scan_ell, ptr_nflag, nsta, ninc, and npix are declared as PRIVATE(static)
'int32_t *longptr' was added in get_l1a_record() routine and used as the
lhs to nav.nflag so that to match the prototype.
 
Revision 1.2  1994/05/10 18:44:37  seawifst
May 6, 1994 version 1.2
 
Revision 1.1  1994/04/19 13:23:38  seawifst
Initial revision
 
Removed all code that was only compiled when GENERAL, USE_STRUCT,
or TESTCODE were defined. (687 lines of code removed.)
Norman Kuring       6-Nov-1996
 
Added code to handle "Calibration" Vgroup data.
Norman Kuring       6-Nov-1996
 
Removed all code that was only compiled when EVERYTHING was defined.
Removed non-prototype style function definitions.
Norman Kuring       7-Nov-1996
 
Put back all of the GENERAL and EVERYTHING code.  Heavy sigh.
Norman Kuring       25-Nov-1996
 
fix values of l1a_data and geoloc when stray light processing is done.
W. Robinson  3-Apr-1997
 
W. Robinson , GSC, 20 Mar 98  update the dark restore calculation to filter
unreasonable values and to take the median
 */
 
 
#include <stdio.h>
#include <string.h>
#include <stdlib.h>
/* for ceil used in STRAY_LIGHT_COEF macro */
#include <math.h>
 
// why reinvent the wheel
#include <gsl/gsl_sort.h>
#include <gsl/gsl_statistics.h>
 
#include "hdf.h"
//#include  "hdfhdr.h"
//#include  "hdfmac.h"
//#include  "usrhdr.h"
#include "usrmac.h"
 
//#include  "generic.h"
 
#include "ffm.h"
#include "cdl_object.h"
#include "hdf_object.h"
#include "SeaWiFS.h"
#include "navigation.h"
#include "l1a.h"
#include "level_1a_index.h"
 
#include "datatype.h"
#include "tlm.h"
#include "l1a_proto.h"
 
/* prototype for get_WIFSinfo() */
#include "WIFSHDF.h"
 
/* stray light */
#include "st_lt.h"
#include "cal_l1a.h"
#include "get_cal.h"
#ifndef NO_SET_HIGHLT1
#include "l2_flags_sw.h"
#endif /* !NO_SET_HIGHLT1 */
 
/* cpp_processing_start_here_please_do_not_delete_this_line */
 
/* Make following as default in all future release(040194)      */
 
 
#define FIRST_GET_L1A_OPEN 1
 
#ifndef DEFAULT_BUFFER_BLKSIZE
#define DEFAULT_BUFFER_BLKSIZE 10
#endif
 
#define BLKSIZE_NOT_CHANGED 0
#define BLKSIZE_CHANGED  1
 
/*
   Following are for get_l1a_record routine to keep track of:
   1. current access record number
   2. total record number
 */
 
#define CLEAR 0
#define INIT 1
#define RECORDS 2
 
#ifndef NO_SET_HIGHLT1
#define FIRST_KNEE 1
/* this global variable is in the code with calibrate_l1a()     */
extern float32 cal_counts[BANDS_DIMS_1A][GAINS_DIMS_1A][KNEES_DIMS_1A];
extern float32 cal_rads[BANDS_DIMS_1A][GAINS_DIMS_1A][KNEES_DIMS_1A];
 
/* global storage space for BAND 8 gain value to be used in stray_light */
short gain8;
 
#endif /* !NO_SET_HIGHLT1 */
 
#define STRAY_LIGHT_COMPLETE 0
#define STRAY_LIGHT_LATE_START 1
#define STRAY_LIGHT_EARLY_STOP 2
#define STRAY_LIGHT_RANDOM 4
#define PIX_SUB_V 4.0
#define GAC_STRAY_LIGHT_COEF(v) (ceil(roundingup = v/PIX_SUB_V))
#define LAC_STRAY_LIGHT_COEF(v) (v)
#define STRAY_LIGHT_COEF(v,dtyp) (strcmp(dtyp,"GAC")?LAC_STRAY_LIGHT_COEF(v):GAC_STRAY_LIGHT_COEF(v))
 
#define NTRY_WARNING 5
#define NTRY_LIMIT 20
 
 
#ifdef L1A_INIT_NEEDED
 
#ifdef PROTOTYPE
PRIVATE void init(void)
#else
 
PRIVATE void init()
#endif
{
    char *FUNC = "l1a_init";
    int i, j;
 
    ON_TIME();
    /* invalidate all fids                      */
    for (i = 0; i < MAX_HDF_L1AGET; i++) {
        fid[i] = -1;
        blksize_changed[i] = BLKSIZE_NOT_CHANGED;
        blksize[i] = DEFAULT_BUFFER_BLKSIZE;
        /* set to zero, indicate no call to get_l1a_record yet  */
        blksrec[i] = 0;
        cal_table_file[i] = NULL;
        stray_light_flag[i] = 0;
        ltyp_frac_coef[i] = 0.0;
        out_band_flag[i] = 0;
        recursive_flag[i] = FALSE;
        recursive_free = 0;
        stray_light_calling_flag[i] = STRAY_LIGHT_COMPLETE;
        stray_light_scan_no[i] = 1;
    }
    OFF_TIME();
}
 
#endif // L1A_INIT_NEEDED
 
/*------------------------------------------------------------------------------
 
    calculate the mean and standard deviation of the 8 bands of 
    dark-restore values for each band
    W. Robinson, GSC, 18 Mar 98  change dark computation to filter
        dark values <5, > 35 and then find the median on the remainder
        Note routines d_stats, i16comp and d_sd are a part of this now
 
------------------------------------------------------------------------------*/
#ifdef PROTOTYPE
int dark_rest_stat(int16 *data, int nrec, float *dark_mean, float *dark_std)
#else
 
int dark_rest_stat(data, nrec, dark_mean, dark_std)
int16 *data; /* 8 bands of dark restore data */
int nrec; /* number of scan lines */
float *dark_mean; /* 8 bands of dark restore mean values */
float *dark_std; /* 8 bands of dark restore std values */
#endif
{
    double *ddata;
    double dmedian;
    double dark;
    int32 num_out, i, j;
    ddata = malloc(nrec * sizeof (double));
 
    for (i = 0; i < 8; i++) {
        for (j = 0; j < nrec; j++) {
            ddata[j] = (double) *(data + j * 8 + i);
        }
        gsl_sort(ddata, 1, nrec);
        dmedian = gsl_stats_median_from_sorted_data(ddata, 1, nrec);
        memset(ddata, 0, nrec);
        num_out = 0;
        for (j = 0; j < nrec; j++) {
            dark = (double) *(data + j * 8 + i);
            if (dark > dmedian * 0.85 && dark < dmedian * 1.15) {
                ddata[num_out] = (double) *(data + j * 8 + i);
                num_out++;
            }
        }
        *(dark_mean + i) = (float) gsl_stats_mean(ddata, 1, num_out);
        *(dark_std + i) = (float) gsl_stats_sd(ddata, 1, num_out);
    }
    free(ddata);
    return 0;
 
}