/*---------------------------------------------------------------------*/
/* fluorescence.c -  functions to support chlorophyll fluorescence.    */
/*                                                                     */
/* It should be noted that these functions are only applicable to a    */
/* sensor with a red band in the range of chlorophyll fluorescence.    */
/* MODIS bands are assumed.                                            */
/*                                                                     */
/* Implementation for MSL12:  B. Franz, Aug 2004.                      */
/*---------------------------------------------------------------------*/

#include "l12_proto.h"

static float pi = PI;

/*---------------------------------------------------------------------*/
/* get_flh - fluorescence line height for each pixel in rec            */
/*           flh output is in radiance units (mW/cm^2/um/sr)           */
/*---------------------------------------------------------------------*/
void get_flh(l2str *l2rec, float flh[])
{
    static float flhmin = 0.0;
    static long  ib667, ib678, ib748;
    static int   firstCall = 1;

    long  ip, ipb;
    float R1;
    float R2;
    float R3;
    float base;

    if (firstCall) {
        firstCall = 0;
        ib667 = windex(667.,l2rec->fwave,NBANDS);
        ib678 = windex(678.,l2rec->fwave,NBANDS);
        ib748 = windex(748.,l2rec->fwave,NBANDS);
    }

    for (ip=0; ip<l2rec->npix; ip++) {

        flh[ip] = 0.0;

        if (!l2rec->mask[ip]) {

            ipb = NBANDS*ip + ib667;
            R1 = l2rec->Lt[ipb] - l2rec->Lr[ipb] - l2rec->TLg[ipb];

            ipb = NBANDS*ip + ib678;
            R2 = l2rec->Lt[ipb] - l2rec->Lr[ipb] - l2rec->TLg[ipb];

            ipb = NBANDS*ip + ib748;
            R3 = l2rec->Lt[ipb] - l2rec->Lr[ipb] - l2rec->TLg[ipb];

            base = R3 + (R1 - R3) * ((746.3 - 676.7) / (746.3 - 665.1));

            flh[ip] = R2 - base;

            if (flh[ip] < flhmin)
	        l2rec->flags[ip] |= CHLRANGE;

	} else {
	    l2rec->flags[ip] |= CHLFAIL;
	}
    }
}


/*---------------------------------------------------------------------*/
/* get_cfe - chlorophyll fluorescence efficiency for each pixel in rec */
/*---------------------------------------------------------------------*/
void get_cfe(l2str *l2rec, float cfe[])
{
    static int   firstCall = 1;
    static float cfemax = 0.1;
    static float *arp;
    static float *flh;

    /* we want flhq in same units as arp (Ein/m^2/s) */
    /* flh is in mW/cm^2/um/sr / 100 = W/m^2/nm/sr    */
    /* hc = 1986.5e-19 J*nm/photon, Lambda 676.7 nm   */
    /* radiance->iradiance = 4*pi                     */
    /* line height -> area = 43.38 nm                 */
    /* 6.023e23 per mole Quanta                       */

    static float flh2flhq = 1./100. * 676.7/(1986.5e-19) * 4.0*PI * 43.38 / 6.023e23;

    float flhq;
    long  ip;

    if (firstCall) {
        firstCall = 0;
        arp = calloc(l2rec->npix, sizeof(float));
        flh = calloc(l2rec->npix, sizeof(float));
    }
    
    get_flh(l2rec,flh);
    get_arp(l2rec,arp);

    for (ip=0; ip<l2rec->npix; ip++) {
        cfe[ip] = 0.0;
        if (!l2rec->mask[ip] && flh[ip] > 0.0 && arp[ip] > 0.0) {
	    flhq = flh[ip]*flh2flhq;
            cfe[ip] = flhq/arp[ip];
	} else {
	    l2rec->flags[ip] |= CHLFAIL;
	}
        if (cfe[ip] > cfemax)
	    l2rec->flags[ip] |= CHLRANGE;
    }
}


/* Notes from modcol (anly8dbl.f90)
! 4*$PI - convert radiance to scalar irradiance
! energy of 1 photon = hc/lambda
! h = 6.6261e-34 [Js]; c = 299,792,458 [m/s] (*10^9 for nm/s)
! energy of one photon = (1986.5*10^-19)/Lambda [W s]
! where lambda is [nm]
!
! FLHQ units are quanta/m^2/s
! to convert FLH into total fluorescence under the fluorescence
! curve :
! The fluorescence curve is described as a gaussian curve with
! center at 683 nm and half maximum emission of 25 nm (from
! Collins et al. 1985
! If L683 (1 nm width) = 1 then the area under the curve= 26.61
! Band14 center= 677 bandwidth= 10 Fluorescence signal read= 8.61
! Band13 center= 665 bandwidth= 10 Fluorescence signal read= 2.86
! Band15 center= 747 bandwidth= 10 Fluorescence signal read= 1.55e-6
! All band readings correspond to a 10 nm bandwidth signal.  For this
! reason the signals must be divided by 10 to get the reading per nm
!
! Because band 13 is affected by the fluorescence signal then there
! is a contribution of this signal equal to 0.2478 to the baseline
! correction.
!
! The conversion factor of FLH into area under the curve can be
! computed: area under the curve / (Band14 read - baseline contrib)
! factor = 26.61 /(0.861 - 0.2478) = 43.38
! rescale: 43.38 (nm) -> 0.04338 (um)
!
! FLHQ = FLH * Lambda/(hc) * (radiance->irradiance)*(FLH->area)
! FLHQ is in quanta m-2 s-1 and ARP is in Einstein m-2 s-1
! 1 Einstein = 1 mol quanta = 6.023e23 quanta
*/