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


#define PI atan2(0.,-1)
#define TPI 2.*atan2(0.,-1)

double a, amp, arg;
double cycles, cy2;
double d, dt, dur;
double f1,f2,freq, fspectral;
double oct;
double rate;
double t, t2;

double chkarg;

unsigned long i;
unsigned long itype;
unsigned long n, ntimes;

char filename[5][25];
FILE *pFile[5];

void main()
{
    printf("\n\n  sineth, ver 1.0 \n");
	printf("\n by Tom Irvine ");
	printf("\n    tomirvine@aol.com ");

	printf("\n This program generates a sine sweep time history.");

	printf("\n");

    strcpy(filename[0],"d.out");
	strcpy(filename[1],"a.out");
	strcpy(filename[2],"f.out");


    pFile[0]=fopen(filename[0],"w");
    pFile[1]=fopen(filename[1],"w");
    pFile[2]=fopen(filename[2],"w");

    printf("\n\n enter starting frequency (Hz)  ");
    scanf("%lf", &f1);

	printf("\n enter ending   frequency (Hz)  ");
    scanf("%lf", &f2);

    printf("\n\n enter amplitude  ");
	scanf("%lf", &amp);

    printf("\n enter duration (sec)  ");
	scanf("%lf",&dur);

	printf("\n enter sweep type:  1=log 2=linear  ");
	scanf("%d",&itype);

    dt=(1./f2)/16.;
	oct=log(f2/f1)/log(2.);

	printf("\n\n number of octaves = %lf \n", oct );
    printf("\n dt = %lf sec",dt);

	ntimes = unsigned long ( dur/dt );

    printf("\n ntimes = %ld ",ntimes );

	cycles=0.;
	t2=dur;

	for(i=0; i<= ntimes; i++)
	{
        t=dt*i;

		if(itype==1)
		{
	
//  weird formula required
		
			if(i==0)
			{
				rate=oct/dur;
				printf("\n\n sweep rate = %lf octaves/sec \n", rate);
			}
			
			fspectral = f1*pow(2.,rate*t);
		
			cy2 = ( f1 /( rate*log(2.) ) )*( -1. + pow(2., (rate*t) ) );

            if(t>0){freq=cy2/t;}

		}
		else
		{

// 0.5 factor is necessary to obtain correct number of cycles.

			 fspectral= (     (f2-f1)*t/dur ) + f1;
		     freq     = ( 0.5*(f2-f1)*t/dur ) + f1;
             cy2=freq*t;
			 
		}
       
		arg=TPI*cy2;

		cycles+=fspectral*dt;

		d=amp*sin(arg);
		a=-( (amp/386.) * ( pow( TPI*fspectral, 2) ) )*sin(arg);

		fprintf(pFile[0],"%lf %lf\n",t, d );
//		fprintf(pFile[1],"%lf %lf\n",t, a );

		fprintf(pFile[2],"%lf %lf %lf \n",t, fspectral, freq);
	}
	printf("\n\n cycles= %lf \n", cycles );
    
	printf("\n Output files: \n");
	printf("\n    d.out - time(sec)  amplitude ");
	printf("\n    f.out - time(sec)  spectral freq (Hz)  sine argument freq (Hz) \n");

	printf("\n Please call the output files into your own plotting program. \n");


	printf("\n Press any key to exit.");
	getch();
}