/**
* Copyright 2005-2007 ECMWF
*
* Licensed under the GNU Lesser General Public License which
* incorporates the terms and conditions of version 3 of the GNU
* General Public License.
* See LICENSE and gpl-3.0.txt for details.
*/

/*
 *
 * Author: Enrico Fucile <enrico.fucile@ecmwf.int>
 *
 *
 */


#include "grib_api_internal.h"

/*
   This is used by make_class.pl

   START_CLASS_DEF
   CLASS      = box
   SUPER      = grib_box_class_gen
   IMPLEMENTS = destroy
   IMPLEMENTS = get_points
   IMPLEMENTS = init
   MEMBERS = double* lats
   MEMBERS = size_t nlats
   MEMBERS = double** lons
   MEMBERS = long* nlons
   MEMBERS = size_t size

   END_CLASS_DEF

 */

/* START_CLASS_IMP */

/*

Don't edit anything between START_CLASS_IMP and END_CLASS_IMP
Instead edit values between START_CLASS_DEF and END_CLASS_DEF
or edit "box.class" and rerun ./make_class.pl

*/


static void init_class              (grib_box_class*);

static int init               (grib_box* box,grib_handle* h,grib_arguments* args);
static grib_points* get_points(grib_box* box, double north, double west, double south,double east,int *err);
static int destroy            (grib_box* box);

typedef struct grib_box_reduced_gaussian{
  grib_box box;
/* Members defined in gen */
/* Members defined in reduced_gaussian */
	double* lats;
	size_t nlats;
	double** lons;
	long* nlons;
	size_t size;
} grib_box_reduced_gaussian;

extern grib_box_class* grib_box_class_gen;

static grib_box_class _grib_box_class_reduced_gaussian = {
    &grib_box_class_gen,                         /* super                     */
    "reduced_gaussian",                         /* name                      */
    sizeof(grib_box_reduced_gaussian),      /* size of instance          */
    0,                              /* inited */
    &init_class,                    /* init_class */
    &init,                          /* constructor               */
    &destroy,                       /* destructor                */
    &get_points,                    /* get points           */
};

grib_box_class* grib_box_class_reduced_gaussian = &_grib_box_class_reduced_gaussian;


static void init_class(grib_box_class* c)
{
}
/* END_CLASS_IMP */

static int init(grib_box* box,grib_handle* h,grib_arguments* args)
{
	grib_box_reduced_gaussian* self=(grib_box_reduced_gaussian*)box;
	int n=1;
	int ret=0;
	char* key;
	long* pl;
	double lat_first,lat_last,lon_first,lon_last;
	long order;
	double *lats;
	long j,l,ilon_first,ilon_last,i;
	long row_count;
	grib_context* c=box->context;

	key = (char*)grib_arguments_get_name(h,args,n++);
	ret=grib_get_double(h,key,&lat_first);  
	if (ret) return ret;

	key = (char*)grib_arguments_get_name(h,args,n++);
	ret=grib_get_double(h,key,&lon_first);  
	if (ret) return ret;

	key = (char*)grib_arguments_get_name(h,args,n++);
	ret=grib_get_double(h,key,&lat_last);  
	if (ret) return ret;

	key = (char*)grib_arguments_get_name(h,args,n++);
	ret=grib_get_double(h,key,&lon_last);  
	if (ret) return ret;

	key = (char*)grib_arguments_get_name(h,args,n++);
	ret=grib_get_long(h,key,&order);  
	if (ret) return ret;

	key = (char*)grib_arguments_get_name(h,args,n++);
	ret=grib_get_size(h,key,&self->nlats);
	if (ret) return ret;
	pl=grib_context_malloc(h->context,self->nlats*sizeof(long));
	ret=grib_get_long_array(h,key,pl,&self->nlats);  
	if (ret) return ret;

	lats=(double*)grib_context_malloc(h->context,sizeof(double)*order*2);
	if((ret = grib_get_gaussian_latitudes(order, lats)) != GRIB_SUCCESS)
	        return ret;

	if (self->nlats==order*2) {
		/*global (latitudes)*/
		self->lats=lats;
	} else {
		/*sub area (latitudes)*/
		int l=0;
		double d=fabs(lats[0]-lats[1]);
		while (fabs(lat_first-lats[l]) > d ) {l++;}
		self->lats=(double*)grib_context_malloc(h->context,sizeof(double)*self->nlats);
		for (i=0;i<self->nlats;i++) self->lats[i]=lats[l++];
		grib_context_free(box->context,lats);
	}

	self->lons=grib_context_malloc_clear(box->context,sizeof(double*)*self->nlats);
	self->size=0;
	if (lon_first != 0                 ||
			fabs(lon_last  - (360.0-90.0/order)) > 90.0/order) {
		/* sub area (longitudes) */
		for (j=0;j<self->nlats;j++) {
			row_count=0;
			grib_get_reduced_row(pl[j],lon_first,lon_last,
                  &row_count,&ilon_first,&ilon_last);
			self->size+=row_count;
			if (ilon_first>ilon_last) ilon_first-=pl[j];
			l=0;
			self->lons[j]=grib_context_malloc_clear(c,sizeof(double)*row_count);
			for (i=ilon_first;i<=ilon_last;i++) 
				self->lons[j][l++]=((i)*360.0)/pl[j];
			pl[j]=row_count;
		}
	} else {
		/* global (longitudes) */
		for (j=0;j<self->nlats;j++) {
			self->lons[j]=grib_context_malloc_clear(c,sizeof(double)*pl[j]);
			self->size+=pl[j];
			for (i=0;i<pl[j];i++) 
				self->lons[j][i]=(i*360.0)/pl[j];
		}

	}
	self->nlons=pl;

	return 0;
}

static int destroy(grib_box* box) {
	grib_box_reduced_gaussian* self=(grib_box_reduced_gaussian*)box;
	grib_context* c=box->context;
	int i;

	grib_points_delete(box->points);

	grib_context_free(c,self->lats);
	for (i=0;i<self->nlats;i++) 
		grib_context_free(c,self->lons[i]);
	grib_context_free(c,self->lons);
	
  return GRIB_SUCCESS;
}

static grib_points* get_points(grib_box* box, double north, double west, double south,double east, int *err) {
	long j,i,index,l;
	double lat,lon;
	size_t group_count;
	int first=1;
	grib_box_reduced_gaussian* self=(grib_box_reduced_gaussian*)box;
	grib_points* points=NULL;

	grib_context* c=box->context;
	points=grib_points_new(c,self->size);
	if (!points) 
		grib_context_log(c,GRIB_LOG_FATAL,"unable to create grib_points\n");

	index=0;
	l=0;
	group_count=0;
	for (j=0;j<self->nlats;j++) {
		lat=self->lats[j];
		for (i=0;i<self->nlons[j];i++)  {
			lon=self->lons[j][i];
			if (lat<north && lat>south && lon>west && lon<east) {
				points->latitudes[l]=lat;
				points->longitudes[l]=lon;
				points->indexes[l]=index;
				points->group_len[group_count]++;
				if (first) {
					points->group_start[group_count]=index;
					points->group_len[group_count]=0;
					first=0;
				} else {
					if (index!=points->group_start[group_count]+points->group_len[group_count]) {
						group_count++;
						points->group_start[group_count]=index;
						points->group_len[group_count]=0;
						points->n_groups=group_count+1;
					} 
				}
				l++;
			}
			index++;
		}
	}
	points->group_len[group_count]++;
	points->n=l;
	if (box->points) grib_points_delete(box->points);
	box->points=points;

	*err=0;

	return points;
}