#include <math.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <time.h>

#define round(x) ((x)>=0?(long)((x)+0.5):(long)((x)-0.5))

/* globale Arbeitsvariablen */

#include "cube.h"

#define FOP 25


double sinc(double t, double T)
/* 
	t time
	T period
*/
{
	if ( t == 0 ) return 1.0;
	else
	{
		return	sin(M_PI*t/T)/(M_PI*t/T);
	}
}
double window(double t, double T, long n)
/*
	Blackmann-Nutall window
	t time
	T period
	n filter length
*/
{
	if (fabs(t)>(n*T)) return 0.0;

	else return 0.3635819+0.4891775*cos(M_PI*t/(n*T))+0.1365995*cos(2*M_PI*t/(n*T))+0.0106411*cos(3*M_PI*t/(n*T));
}



void main(argc,argv)
int argc;
char **argv;
{
  long			nsamp1,nsamp2,j,ch,jj;
  long			tqx;
  double		tx,tx2,t1,t2,t,T; 
  long			ts1,tmks1,ts2,tmks2,hx,gps_t,trs,trmks;
  long			delta_t,tt;

  char  		buf[1000];
  int			code;
  int			ampx[FOP+FOP+1][10];
  double		yx1[10];
  long long		ts1l[FOP+FOP+1],tsrl;
  double		coeff[(FOP+1)*40000+1];


extern int optind, opterr, optopt;
   int c,errflg; 
   time_t start_time;
   struct tm *tss;

   errflg=0;
   putenv ("TZ=GMT");
   tzset; 

  gps_t=0;
  opterr=0; while ((c = getopt (argc, argv, "o:t")) != -1)
  switch (c)
      {
         case 't':        /* gps_t timing info explicite */
           gps_t=1;               
           break;
         case '?':
           fprintf(stderr,"ERROR: Unrecognized option: -%c\n", argv[0], optopt);
           errflg++;
      }

/* read config */   
   if (fgets(buf,1000,stdin)) 
   {
	code=sscanf(buf,STR_FORMAT,(char *) SC_STR);
   }

   fprintf(stdout,STR_FORMAT"\n",PR_STR);


/* read first sample */
for (j=0; j<(FOP+FOP); j++)
{
  if (fgets(buf,1000,stdin))
  {
	switch (CH_NUM)
	{
		case 1: code=sscanf(buf,"%ld%d%ld%ld",&nsamp1,&ampx[j][0],&ts1,&tmks1); break;
		case 2: code=sscanf(buf,"%ld%d%d%ld%ld",&nsamp1,&ampx[j][0],&ampx[j][1],&ts1,&tmks1); break;
		case 3: code=sscanf(buf,"%ld%d%d%d%ld%ld",&nsamp1,&ampx[j][0],&ampx[j][1],&ampx[j][2],&ts1,&tmks1); break;
	}
	ts1l[j]=ts1*1000000+tmks1;
/*        fprintf(stderr,"get sample %ld time %ld \n", j, ts1l[j]);	*/
/* no data treatment */
	if (nsamp1==-1)
	{ 
		fprintf(stdout,"%ld ",-1L); 
  		for (ch=0; ch<CH_NUM; ch++) fprintf(stdout,"%d ",0); 
  		fprintf(stdout,"%ld %ld\n",0L,0L); 
		fprintf(stderr,"%s ERROR: No useful data in resamp!\n",argv[0]); 
		exit(77);
	}
  }
}



/* sample increment in mks */
  hx=(long) floor((1000000.0/S_RATE));

/* construct sinc-window table */
  if (hx>20000)
  { 
  	fprintf(stdout,"%ld ",-1L); 
  	for (ch=0; ch<CH_NUM; ch++) fprintf(stdout,"%d ",0); 
   	fprintf(stdout,"%ld %ld\n",0L,0L); 
 	fprintf(stderr,"%s sample rate of %ld not supported in resamp!\n",argv[0],1000000/hx); 
 	exit(77);
  }

  for (j=0; j<=((FOP+1)*hx*2); j++) 
  {
	coeff[j]=sinc(((double) j), (double) hx)*window(((double) j),(double) hx,FOP);
  }
/* fprintf(stderr,"XXX %ld %ld %ld \n",hx, (FOP+1)*40000+1,((FOP+1)*hx*2)); */



/* start value for new time scale */
  trs=(long) (ts1l[FOP]/1000000);
  tmks1=(long) (ts1l[FOP]-trs*1000000);
  trmks=(long) (ceil((double) (tmks1/(double) hx))*((double)hx));
  if (trmks==1000000) {trmks=0; trs=trs+1;}
  tsrl=trs*1000000+trmks;


/* fprintf(stdout,"Start Value %ld %6ld %ld\n",trs,trmks,ts1l[FOP]); */


/* main loop */
  while (!feof(stdin))
  {

  	if (fgets(buf,1000,stdin))
  	{
		switch (CH_NUM)
		{
			case 1: code=sscanf(buf,"%ld%d%ld%ld",&nsamp2,&ampx[FOP+FOP][0],&ts2,&tmks2); break;		
			case 2: code=sscanf(buf,"%ld%d%d%ld%ld",&nsamp2,&ampx[FOP+FOP][0],&ampx[FOP+FOP][1],&ts2,&tmks2); break;
			case 3: code=sscanf(buf,"%ld%d%d%d%ld%ld",&nsamp2,&ampx[FOP+FOP][0],&ampx[FOP+FOP][1],&ampx[FOP+FOP][2],&ts2,&tmks2); break;
		}
  	}

  	ts1l[FOP+FOP]=ts2*1000000+tmks2;

/* fprintf(stdout,"XXX %ld %ld %ld %ld --- %ld %ld %ld --- ",(long)((ts1l[0]-tsrl)/hx),(long)((ts1l[FOP+FOP]-tsrl)/hx), ts1l[0], ts1l[FOP+FOP], ts1l[FOP],tsrl,ts1l[FOP+1]); */


/* remember: lower time boundary included !!!*/

/* check if new sample is between both (excl. upper) samples */
	while (tsrl<ts1l[FOP+1])
	{
/* sinc interpolation */
		delta_t=ts1l[FOP]-tsrl;

		for (ch=0; ch<CH_NUM; ch++)
		{
			yx1[ch]=0.0;
			for (jj=0; jj<=(2*FOP) ; jj++)
			{
				tt=abs((jj-FOP)*hx+delta_t); 
/*				tt=abs(ts1l[jj]-ts1l[FOP]+delta_t); */

/*				yx1[ch]=yx1[ch]+sinc(tt, hx)*window(tt,hx,FOP)*ampx[jj][ch]; */
				yx1[ch]=yx1[ch]+coeff[tt]*(double) (ampx[jj][ch]); 
			}
		} 
		if (gps_t==0) 
		{
			fprintf(stdout,"%ld %6ld ",(long) (tsrl/1000000),(long)(tsrl - ((long) (tsrl/1000000))*1000000));
		}
		else 
		{
			start_time=(long) (tsrl/1000000);
  			tss = localtime(&start_time); 
			fprintf(stdout,"%4d-%02d-%02dT%02d:%02d:%02d %6ld ",(int) (tss->tm_year)+1900,(int) (tss->tm_mon)+1,(int) (tss->tm_mday),(int) (tss->tm_hour),(int) (tss->tm_min),(int) (tss->tm_sec),(long)(tsrl - ((long) (tsrl/1000000))*1000000)); 
		}

/* TEST sin 
		t=tsrl*0.000001-1400054271.0;
		T=10.0;
		for (ch=0; ch<CH_NUM; ch++) yx1[ch] = (int) round((double) yx1[ch])-(int) (1000000.0*sin(M_PI*t/T));
*/

		for (ch=0; ch<CH_NUM; ch++) fprintf(stdout," %d",(int) round((double) yx1[ch])); fprintf(stdout,"\n");


/* update new time */
		tsrl=tsrl+hx; 
	}


/* update data vector */
	for (j=1; j<=(FOP+FOP); j++) 
	{
		ts1l[j-1]=ts1l[j];
		for (ch=0; ch<CH_NUM; ch++) ampx[j-1][ch]=ampx[j][ch];
	}


  }
 }


