41     if( ( 
l1rec->cld_dat = malloc( 
sizeof( cld_struc ) ) ) == 
NULL )
 
   43       printf( 
"-E- %s, %d: Failure to allocate cloud albedo data structure\n",
 
   47     if( ( ( 
l1rec->cld_dat->sfc_albedo_659 = (
float *) 
 
   48             malloc( 
npix * 
sizeof( 
float ) ) ) == 
NULL ) || 
 
   49         ( ( 
l1rec->cld_dat->sfc_albedo_858 = (
float *)
 
   50             malloc( 
npix * 
sizeof( 
float ) ) ) == 
NULL ) ||
 
   51         ( ( 
l1rec->cld_dat->sfc_albedo_1240 = (
float *)
 
   52             malloc( 
npix * 
sizeof( 
float ) ) ) == 
NULL ) ||
 
   53         ( ( 
l1rec->cld_dat->sfc_albedo_1640 = (
float *)
 
   54             malloc( 
npix * 
sizeof( 
float ) ) ) == 
NULL ) ||
 
   55         ( ( 
l1rec->cld_dat->sfc_albedo_2130 = (
float *)
 
   56             malloc( 
npix * 
sizeof( 
float ) ) ) == 
NULL ) ||
 
   57         ( ( 
l1rec->cld_dat->cth_alb_init = (
float *)
 
   58             malloc( 
npix * 
sizeof( 
float ) ) ) == 
NULL ) ||
 
   59         ( ( 
l1rec->cld_dat->cth_alb_unc_init = (
float *)
 
   60             malloc( 
npix * 
sizeof( 
float ) ) ) == 
NULL ) )
 
   62       printf( 
"-E- %s, %d: Failure to allocate cloud albedo data structures\n",
 
   69 int read_albedo( int32_t d_np, int32_t d_nl, int32_t st_lin, int32_t ix_clim, 
 
   70   int ncid, 
int *d_id, 
unsigned char ***alb_dat )
 
   98   int32_t iwav, i_np, i_nl;
 
  101   unsigned char *xfr_arr, **alb_lcl;
 
  105   if( alb_lcl[0] != 
NULL )
 
  107     for( iwav = 0; iwav < 5; iwav++)
 
  108       free( alb_lcl[iwav] );
 
  111   if( ( xfr_arr =  (
unsigned char *) malloc( d_np * d_nl * 
 
  112      sizeof( 
unsigned char ) ) ) == 
NULL )
 
  114      printf( 
"-E- %s, %d: Error allocating albedo read storage\n", 
 
  115        __FILE__, __LINE__ );
 
  123   for( iwav = 0; iwav < 5; iwav++)
 
  125     if( ( alb_lcl[iwav] = (
unsigned char *)
 
  126       malloc( i_np * i_nl * 
sizeof( 
unsigned char ) ) ) == 
NULL )
 
  128       printf( 
"-E- %s, %d: Error allocating albedo read storage\n", 
 
  129         __FILE__, __LINE__ );
 
  151     if( nc_get_vara_uchar( ncid, d_id[iwav], 
start, 
count, xfr_arr ) != 
 
  154       printf( 
"-E- %s, %d: Error reading albedo\n",
 
  155        __FILE__, __LINE__ );
 
  160     for( ilat = 0; ilat < d_nl; ilat++ )
 
  162       memcpy( ( alb_lcl[iwav] + ilat * i_np + 1 ),
 
  163         ( xfr_arr + ilat * d_np ), d_np );
 
  165       *( alb_lcl[iwav] + ilat * i_np ) =
 
  166         *( xfr_arr + ilat * d_np + ( d_np - 1 ) );
 
  169     memcpy( ( alb_lcl[iwav] + ( i_nl - 1 ) * i_np ), 
 
  170       ( alb_lcl[iwav] + ( i_nl - 2 ) * i_np ), i_np );
 
  201   static int32_t firstcall = 0, is_clim;
 
  202   static int32_t subset_read = 0, scan_npix;
 
  203   static float bdy_siz = 5.;
 
  204   static float scale[5];
 
  205   static float lat_min, lon_min, res1, res2;
 
  206   static double lat_res, lon_res, alb_st_lat, alb_lat_rng[2] = { 900., -900. };
 
  207   static unsigned char **alb_dat, fillv[5];
 
  208   static int32_t ix_clim = -1, sub_nl;
 
  209   static size_t d_np, d_nl, d_nt, i_np;
 
  210   int32_t iwav, ipx, sub_st, good_nav, ix_lat, ix_lon;
 
  211   int32_t n_done, ip, il;
 
  212   float dat_sub[2][2], fin_val, sum;
 
  213   float lat_rng[2] = { 900., -900. };
 
  215   double sec, frac_lat, frac_lon, 
lat, 
lon;
 
  217   char tit_clim[] = 
"MODIS/TERRA+Aqua BRDF/Albedo Gap-Filled Snow-Free 8-Day climatology";
 
  218   char tit_daily[] = 
"MODIS/TERRA+Aqua BRDF/Albedo Gap-Filled Snow-Free data";
 
  219   char str_tit[FILENAME_MAX];
 
  220   static int ncid, d_id[5], dim_id, dim_id2;
 
  223   char *alb_wav_names[5] = { 
"Albedo_Map_0.659", 
"Albedo_Map_0.858", 
 
  224     "Albedo_Map_1.24", 
"Albedo_Map_1.64", 
"Albedo_Map_2.13" };
 
  230     alb_file = 
input->sfc_albedo;
 
  232     if (nc_open(alb_file, 0, &ncid) != NC_NOERR) {
 
  234         "-E- %s %d: file: %s is not netcdf, not acceptable albedo file\n",
 
  235                     __FILE__, __LINE__, alb_file);
 
  239     nc_get_att_text(ncid, NC_GLOBAL, 
"title", str_tit);
 
  240     if( strncmp( str_tit, tit_daily, strlen(tit_daily) ) == 0 )
 
  242       printf( 
"%s, %d: TEMP: we found the daily file: %s\n", 
 
  243         __FILE__, __LINE__, str_tit );
 
  246     else if( strncmp( str_tit, tit_clim, strlen(tit_clim) ) == 0 )
 
  248       printf( 
"%s, %d: TEMP: we found the clim file: %s\n", __FILE__, 
 
  249         __LINE__, tit_clim );
 
  253       ix_clim = ( doy - 1 ) / 8;  
 
  254       printf( 
"-I- %s, %d: Using climate doy # %d, index: %d\n", __FILE__, 
 
  255         __LINE__, doy, ix_clim );
 
  256       printf( 
"-I- From file: %s\n", alb_file );
 
  260       printf( 
"-E- %s, %d: Incorrect title in albedo dataset: %s\n",
 
  261         __FILE__, __LINE__, alb_file);
 
  266     if( ( nc_get_att_float(ncid, NC_GLOBAL, 
"geospatial_lat_min", &lat_min ) 
 
  268         ( nc_get_att_float(ncid, NC_GLOBAL, 
"geospatial_lon_min", &lon_min )  
 
  271       printf( 
"-E- %s, %d: Unable to read geolocation information for albedo file: %s\n",
 
  272         __FILE__, __LINE__, alb_file);
 
  278     nc_inq_atttype( ncid, NC_GLOBAL, 
"geospatial_lat_resolution", &att_typ );
 
  279     if( att_typ == NC_FLOAT ) {
 
  280       nc1 = nc_get_att_float(ncid, NC_GLOBAL, 
"geospatial_lat_resolution",
 
  283       nc1 = nc_get_att(ncid, NC_GLOBAL, 
"geospatial_lat_resolution",
 
  289     nc_inq_atttype( ncid, NC_GLOBAL, 
"geospatial_lon_resolution", &att_typ );
 
  290     if( att_typ == NC_FLOAT ) {
 
  291       nc2 = nc_get_att_float(ncid, NC_GLOBAL, 
"geospatial_lon_resolution",
 
  294       nc1 = nc_get_att(ncid, NC_GLOBAL, 
"geospatial_lon_resolution",
 
  299     if( ( nc1 != NC_NOERR ) || ( nc2 != NC_NOERR ) )
 
  301       printf( 
"-E- %s, %d: Unable to read geolocation information for albedo file: %s\n",
 
  302         __FILE__, __LINE__, alb_file);
 
  309     scan_npix = 
l1rec->npix;
 
  311     if ((nc_inq_dimid(ncid, 
"lat", &dim_id) != NC_NOERR) ||
 
  312         (nc_inq_dimid(ncid, 
"lon", &dim_id2) != NC_NOERR) )
 
  314       printf( 
"-E- %s, %d: Error retrieving the lat/lon size\n", 
 
  315         __FILE__, __LINE__ );
 
  318     if ( ( nc_inq_dimlen(ncid, dim_id, &d_nl ) != NC_NOERR) ||
 
  319          ( nc_inq_dimlen(ncid, dim_id2, &d_np ) != NC_NOERR) )
 
  321       printf( 
"-E- %s, %d: Error retrieving the lat/lon size\n", 
 
  322         __FILE__, __LINE__ );
 
  327       if (nc_inq_dimid(ncid, 
"doy", &dim_id) != NC_NOERR)
 
  329         printf( 
"-E- %s, %d: Error retrieving the doy size\n", 
 
  330           __FILE__, __LINE__ );
 
  333       if ( nc_inq_dimlen(ncid, dim_id, &d_nt ) != NC_NOERR)
 
  335         printf( 
"-E- %s, %d: Error retrieving the doy size\n", 
 
  336           __FILE__, __LINE__ );
 
  342      alb_dat = (
unsigned char **) malloc( 5 * 
sizeof( 
unsigned char * ) );
 
  348     if( (
long)d_np * (long)d_nl < 4000000 ) subset_read = 0;
 
  351     for( iwav = 0; iwav < 5; iwav++ )
 
  355       if( nc_inq_varid(ncid, alb_wav_names[iwav], ( d_id + iwav ) ) != NC_NOERR)
 
  357         printf( 
"-E- %s, %d: Error retrieving the band id\n", 
 
  358           __FILE__, __LINE__ );
 
  362       if( ( nc_get_att_uchar( ncid, d_id[iwav], 
"_FillValue", 
 
  363         ( fillv + iwav ) ) != NC_NOERR ) ||
 
  364           ( nc_get_att_float( ncid, d_id[iwav], 
"scale_factor",
 
  365         ( 
scale + iwav ) ) != NC_NOERR ) )
 
  367         printf( 
"-E- %s, %d: Error retrieving the fill or scale  value\n", 
 
  368          __FILE__, __LINE__ );
 
  372     if( subset_read == 0 )
 
  374       if( 
read_albedo( d_np, d_nl, 0, ix_clim, ncid, d_id, &alb_dat ) != 0 )
 
  392   if( subset_read == 1 ) 
 
  395     for( ipx = 0; ipx < scan_npix; ipx++ )
 
  398       if( ( 
lat >= -90. ) && ( 
lat <= 90. ) )
 
  400         if( 
lat < lat_rng[0] ) lat_rng[0] = 
lat;
 
  401         if( 
lat > lat_rng[1] ) lat_rng[1] = 
lat;
 
  405     if( ( alb_lat_rng[0] > lat_rng[0] ) || ( alb_lat_rng[1] < lat_rng[1] ) )
 
  408       alb_lat_rng[0] = lat_rng[0] - bdy_siz;
 
  409       if( alb_lat_rng[0] < -90. ) 
 
  410         alb_lat_rng[0] = -90.;
 
  413         sub_st = ( alb_lat_rng[0] + 90. ) / lat_res;
 
  414         alb_lat_rng[0] = ( sub_st * lat_res ) - 90.;
 
  416       alb_lat_rng[1] = lat_rng[1] + bdy_siz;
 
  417       if( alb_lat_rng[1] > 90. ) 
 
  418         alb_lat_rng[1] = 90.;
 
  420         sub_st = ( alb_lat_rng[1] + 90. ) / lat_res;
 
  421         alb_lat_rng[1] = ( sub_st * lat_res ) - 90.;
 
  423       alb_st_lat = alb_lat_rng[0];
 
  427       printf( 
"Latitude: %f - %f\n\n\n",
 
  428         alb_lat_rng[0], alb_lat_rng[1] );
 
  429       sub_nl = ( alb_lat_rng[1] - alb_lat_rng[0] + lat_res/2. ) / lat_res;
 
  430       sub_st = ( alb_st_lat + 90. ) / lat_res;
 
  431       if( 
read_albedo( d_np, sub_nl, sub_st, ix_clim, ncid, d_id, 
 
  437   for( ipx = 0; ipx < scan_npix; ipx++ )
 
  444     if( ( 
lat >= -90. ) && ( 
lat <= 90. ) 
 
  445       && ( 
lon >= -180. ) && (
lon <= 180. ) )
 
  448       frac_lat = ( 
lat - alb_st_lat ) / lat_res;
 
  449       frac_lon = ( 
lon + 180 ) / lon_res;
 
  450       ix_lat = (int32_t) frac_lat;
 
  451       ix_lon = (int32_t) frac_lon;
 
  452       frac_lat = frac_lat - (
float) ix_lat; 
 
  453       frac_lon = frac_lon - (
float) ix_lon;
 
  456     for( iwav = 0; iwav < 5; iwav++ )
 
  463         for( ip = 0; ip < 2; ip++ )
 
  465           for( il = 0; il < 2; il++ )
 
  467             if( ix_lat + il >= sub_nl )
 
  468               dat_sub[ip][il] = fillv[iwav];
 
  470               dat_sub[ip][il] = *( alb_dat[iwav] + ( ix_lon + ip ) +
 
  471                 i_np * ( ix_lat + il ) );
 
  472             if( dat_sub[ip][il] != fillv[iwav] )
 
  474               sum += dat_sub[ip][il] * 
scale[iwav];
 
  486           for( ip = 0; ip < 2; ip++ )
 
  487             for( il = 0; il < 2; il++ )
 
  489               if( dat_sub[ip][il] == fillv[iwav] )
 
  490                 dat_sub[ip][il] = sum;
 
  492                 dat_sub[ip][il] *= 
scale[iwav];
 
  495           fin_val = ( 1 - frac_lon ) * ( ( dat_sub[0][0] * ( 1 - frac_lat ) ) +
 
  496                                      ( dat_sub[0][1] * frac_lat ) ) +
 
  497                           frac_lon *   ( ( dat_sub[1][0] * ( 1 - frac_lat ) ) +
 
  498                                      ( dat_sub[1][1] * frac_lat  ) );
 
  504         case 0: 
l1rec->cld_dat->sfc_albedo_659[ipx] = fin_val;
 
  506         case 1: 
l1rec->cld_dat->sfc_albedo_858[ipx] = fin_val;
 
  508         case 2: 
l1rec->cld_dat->sfc_albedo_1240[ipx] = fin_val;
 
  510         case 3: 
l1rec->cld_dat->sfc_albedo_1640[ipx] = fin_val;
 
  512         case 4: 
l1rec->cld_dat->sfc_albedo_2130[ipx] = fin_val;
 
  543   static size_t start[] = { 0,0,0,0}, 
count[] = { 0,0,0,0};
 
  544   static float *alb_arr, *alb_unc_arr;
 
  545   static int32_t firstcall = 0;
 
  546   float *arr_ptr, talb;
 
  547   size_t nwave, nmon, nlat, nlon, ix_clim;
 
  548   int16_t year, mon, 
day, hh, mm;
 
  549   int32_t ipix, ix, iy;
 
  551   char *alb_file, str_tit[5000];
 
  552   char tit_clim[] = 
"Surface Lambertian-equivalent reflectivity (LER) observed by TROPOMI";
 
  553   char *arr_nms[] = { 
"minimum_LER_clear", 
"uncertainty_clear" };
 
  554   int ncid, d_id, dim_id[4];
 
  560     alb_file = 
input->cth_albedo;
 
  562     printf( 
"%s, %d: I - Loading Cloud Top Height albedo file: %s\n", 
 
  563       __FILE__, __LINE__, alb_file );
 
  564     if (nc_open(alb_file, 0, &ncid) != NC_NOERR) {
 
  566         "-E- %s %d: file: %s is not netcdf, not acceptable cth albedo file\n",
 
  567                     __FILE__, __LINE__, alb_file);
 
  573     status = nc_get_att_text(ncid, NC_GLOBAL, 
"title", str_tit);
 
  583     if( strncmp( str_tit, tit_clim, strlen(tit_clim) ) == 0 )
 
  585       printf( 
"%s, %d: TEMP: we found the cth albedo clim file: %s\n", 
 
  586         __FILE__, __LINE__, tit_clim );
 
  590       printf( 
"-I- %s, %d: Using climate month # %d, index: %ld\n", __FILE__, 
 
  591         __LINE__, mon, ix_clim );
 
  592       printf( 
"-I- From file: %s\n", alb_file );
 
  596       printf( 
"-E- %s, %d: Incorrect title in albedo dataset: %s\n",
 
  597         __FILE__, __LINE__, alb_file);
 
  602     if( ( nc_inq_dimid( ncid, 
"wavelength", dim_id ) != NC_NOERR) ||
 
  603         ( nc_inq_dimid( ncid, 
"month", ( dim_id + 1 ) ) != NC_NOERR) ||
 
  604         ( nc_inq_dimid( ncid, 
"latitude", ( dim_id + 2 ) ) != NC_NOERR) ||
 
  605         ( nc_inq_dimid( ncid, 
"longitude", ( dim_id + 3 ) ) != NC_NOERR) )
 
  607       printf( 
"-E- %s, %d: Error retrieving the wave, month, lat, lon size\n", 
 
  608         __FILE__, __LINE__ );
 
  611     if ( ( nc_inq_dimlen(ncid, dim_id[0], &nwave ) != NC_NOERR) ||
 
  612          ( nc_inq_dimlen(ncid, dim_id[1], &nmon ) != NC_NOERR) ||
 
  613          ( nc_inq_dimlen(ncid, dim_id[2], &nlat ) != NC_NOERR) ||
 
  614          ( nc_inq_dimlen(ncid, dim_id[3], &nlon ) != NC_NOERR) )
 
  616       printf( 
"-E- %s, %d: Error retrieving the wave, month, lat, lon size\n", 
 
  617         __FILE__, __LINE__ );
 
  624       printf( 
"-E- %s, %d: cth albedo file dimension sizes, wave, month, lat, lon are not as expected\n",
 
  625         __FILE__, __LINE__ );
 
  638     if( ( ( alb_arr = (
float *) malloc( nlat * nlon * 
sizeof( 
float ) ) ) 
 
  640         ( ( alb_unc_arr = (
float *) malloc( nlat * nlon * 
sizeof( 
float ) ) )
 
  643       printf( 
"-E- %s, %d: Unable to allocate cth albedo storage\n", 
 
  644         __FILE__, __LINE__ );
 
  648     for( ipix = 0; ipix < 2; ipix++ )
 
  650       if( nc_inq_varid( ncid, arr_nms[ipix], &d_id ) != NC_NOERR )
 
  652         printf( 
"-E- %s, %d: Unable to find dataset %s in cth albedo file\n",
 
  653           __FILE__, __LINE__, arr_nms[ipix] );
 
  656       arr_ptr = ( ipix == 0 ) ? alb_arr : alb_unc_arr;
 
  660         printf( 
"-E- %s, %d: Unable to read dataset %s in cth albedo file\n",
 
  661           __FILE__, __LINE__, arr_nms[0] );
 
  670   for( ipix = 0; ipix < 
l1rec->npix; ipix++ )
 
  672     iy = ( 
l1rec->lat[ipix] + 90. ) / 0.125;
 
  673     ix = ( 
l1rec->lon[ipix] + 180. ) / 0.125;
 
  674     iy = ( iy < 0 )? 0 : iy;
 
  675     iy = ( iy > 1439 ) ? 1439 : iy;
 
  676     ix = ( ix < 0 ) ? 0 : ix;
 
  677     ix = ( ix > 2879 ) ? 2879 : ix;
 
  679     talb = *( alb_arr + iy + 1440 * ix );
 
  680     l1rec->cld_dat->cth_alb_init[ipix] = ( talb > 0.99 ) ? 0.99 : talb;
 
  681     talb = *( alb_unc_arr + iy + 1440 * ix );
 
  682     l1rec->cld_dat->cth_alb_unc_init[ipix] = ( talb < 0.02 ) ? 0.02 : talb;