#include "../../public/ccl/inc/gendefs.h"

#define NBR_PTS 32
#define SL_ZASTINENI_ID 1001
#define SL_APODIZACE_ID 1002

cmd
zastineni()
{
/* V promenne "zastin" je ulozen polomer centralniho zastineni */
 
double souradnice_x1[2*NBR_PTS + 1];
double souradnice_y1[2*NBR_PTS + 1];
double souradnice_x2[2*NBR_PTS + 1];
double souradnice_y2[2*NBR_PTS + 1];

double zastin;
 zastin=0.00001;

double x; //souradnice v obrazove rovine
double pom1,pom2,pom3;
int i;

/*Vypocet rozdeleni intenzity obrazu bodu pro kruhovou pupilu v definovanych bodech */ 
	for (i = -NBR_PTS; i <= NBR_PTS; i++) //vypocet v 2*NBR_PTS bodech
	{
	  x = 4*pi*i/NBR_PTS; //zakresleni probehne v rozsahu -4*pi az 4*pi
	  souradnice_x1[i + NBR_PTS] = x;  //souradnice x se ulozi do pole
 
	  if(x!=0)
	    souradnice_y1[i + NBR_PTS] = (4*j1(x)*j1(x))/(x*x); //vypoctene hodnoty intenzity se ulozi do pole souradnice_y1[]
	  else
	    souradnice_y1[i+NBR_PTS]=1;
	}
	
/*Vypocet rozdeleni intenzity obrazu bodu pro pupilu ve tvaru mezikruzi v definovanych bodech */ 
	for (i = -NBR_PTS; i <= NBR_PTS; i++)
	{
	  x = 4*pi*i/NBR_PTS; //zakresleni probehne v rozsahu -4*pi az 4*pi
	  souradnice_x2[i + NBR_PTS] = x; //souradnice x se ulozi do pole

	  pom1=j1(x);
	  pom2=j1(x*zastin);

	  if(x!=0)
	    { 
	      pom3= (2*pom1/x - (zastin*zastin)*2*pom2/(x*zastin));
	      souradnice_y2[i + NBR_PTS] = pom3*pom3/((1-zastin*zastin)*(1-zastin*zastin)) ; //vypoctene hodnoty intenzity se ulozi do pole souradnice_y2[]
	    }
	      else 
		souradnice_y2[i+NBR_PTS]=1;
	    
	}

	
	// nasleduje kresleni
	graphwin_open(SLIDERBASEWIN);
	graphwin_reset();
graphwin_sliderassign(SL_ZASTINENI_ID, real, drag, "zastineni", zastineni_CB,0,0.9999,0);
	graphwin_slidershow();

	gwt "Rozdeleni intenzity obrazu bodu";
	gvp(0.0, 1.0, 0.0, 1.0, 1.0);
	epset;

	gplot("2", 2*NBR_PTS + 1,
	      souradnice_x1 , souradnice_y1 , "kruhova pupila",
	      souradnice_x2, souradnice_y2," mezikruzi",
	      "relativni intenzita", "x=pi*y'/(lambda*c)", y, y, scs, 0.0, 0.0, -4*pi, +4*pi);
}

cmd zastineni_CB(int slider_id, double zastin)
{
/* V promenne "zastin" je ulozen polomer centralniho zastineni */
  if (zastin==0) zastin=0.0001;

double souradnice_x1[2*NBR_PTS + 1];
double souradnice_y1[2*NBR_PTS + 1];
double souradnice_x2[2*NBR_PTS + 1];
double souradnice_y2[2*NBR_PTS + 1];

double x; //souradnice v obrazove rovine
double pom1,pom2,pom3;
int i;

/*Vypocet rozdeleni intenzity obrazu bodu pro kruhovou pupilu v definovanych bodech */ 
	for (i = -NBR_PTS; i <= NBR_PTS; i++) //vypocet v 2*NBR_PTS bodech
	{
	  x = 4*pi*i/NBR_PTS; //zakresleni probehne v rozsahu -4*pi az 4*pi
	  souradnice_x1[i + NBR_PTS] = x;  //souradnice x se ulozi do pole
 
	  if(x!=0)
	    souradnice_y1[i + NBR_PTS] = (4*j1(x)*j1(x))/(x*x); //vypoctene hodnoty intenzity se ulozi do pole souradnice_y1[]
	  else
	    souradnice_y1[i+NBR_PTS]=1;
	}
	
/*Vypocet rozdeleni intenzity obrazu bodu pro pupilu ve tvaru mezikruzi v definovanych bodech */ 
	for (i = -NBR_PTS; i <= NBR_PTS; i++)
	{
	  x = 4*pi*i/NBR_PTS; //zakresleni probehne v rozsahu -4*pi az 4*pi
	  souradnice_x2[i + NBR_PTS] = x; //souradnice x se ulozi do pole

	  pom1=j1(x);
	  pom2=j1(x*zastin);

	  if(x!=0)
	    { 
	      pom3= (2*pom1/x - (zastin*zastin)*2*pom2/(x*zastin));
	      souradnice_y2[i + NBR_PTS] = pom3*pom3/((1-zastin*zastin)*(1-zastin*zastin)) ; //vypoctene hodnoty intenzity se ulozi do pole souradnice_y2[]
	    }
	      else 
		souradnice_y2[i+NBR_PTS]=1;
	    
	}
	
	// nasleduje kresleni
	graphwin_open(SLIDERBASEWIN);
	eperase;
	gplot("2", 2*NBR_PTS + 1,
	      souradnice_x1 , souradnice_y1 , "kruhova pupila",
	      souradnice_x2, souradnice_y2," mezikruzi",
	      "relativni intenzita", "x=pi*y'/(lambda*c)", y, y, scs, 0.0, 0.0, -4*pi, +4*pi);
}


cmd
apodizace()
{
/* V promenne "apodiz" je ulozen koeficient urcujici prubeh zmeny propustnosti pupily */
 
double souradnice_x1[2*NBR_PTS + 1];
double souradnice_y1[2*NBR_PTS + 1];
double souradnice_x2[2*NBR_PTS + 1];
double souradnice_y2[2*NBR_PTS + 1];

double apodiz;
 apodiz=0.0;

double x; //souradnice v obrazove rovine
double pom1,pom2,pom3;
int i;

/*Vypocet rozdeleni intenzity obrazu bodu pro kruhovou pupilu v definovanych bodech */ 
	for (i = -NBR_PTS; i <= NBR_PTS; i++) //vypocet v 2*NBR_PTS bodech
	{
	  x = 4*pi*i/NBR_PTS; //zakresleni probehne v rozsahu -4*pi az 4*pi
	  souradnice_x1[i + NBR_PTS] = x;  //souradnice x se ulozi do pole
 
	  if(x!=0)
	    souradnice_y1[i + NBR_PTS] = (4*j1(x)*j1(x))/(x*x); //vypoctene hodnoty intenzity se ulozi do pole souradnice_y1[]
	  else
	    souradnice_y1[i+NBR_PTS]=1;
	}
	
/*Vypocet rozdeleni intenzity obrazu bodu pro pupilu jejiz propustnost je promenna dle vztahu g(q)= 1-konst*q    */ 
	for (i = -NBR_PTS; i <= NBR_PTS; i++)
	{
	  x = 4*pi*i/NBR_PTS; //zakresleni probehne v rozsahu -4*pi az 4*pi
	  souradnice_x2[i + NBR_PTS] = x; //souradnice x se ulozi do pole

	  pom1=j1(x);
	  pom2=j0(x);
	 
	  if(x!=0)
	    { 
	      pom3=(2.0/(x*x))*((x*(1.0-apodiz)+4*apodiz/x)*pom1-2*apodiz*pom2);
	      souradnice_y2[i + NBR_PTS] = pom3*pom3; //vypoctene hodnoty intenzity se ulozi do pole souradnice_y2[]
	    }
	      else 
		souradnice_y2[i+NBR_PTS]=1;
	    
	}

	
	// nasleduje kresleni
	graphwin_open(SLIDERBASEWIN);
	graphwin_reset();
graphwin_sliderassign(SL_apodizace_ID, real, drag, "apodizace", apodizace_CB,0,0.9999,0);
	graphwin_slidershow();

	gwt "Rozdeleni intenzity obrazu bodu";
	gvp(0.0, 1.0, 0.0, 1.0, 1.0);
	epset;

	gplot("2", 2*NBR_PTS + 1,
	      souradnice_x1 , souradnice_y1 , "kruhova pupila",
	      souradnice_x2, souradnice_y2," mezikruzi",
	      "relativni intenzita", "x=pi*y'/(lambda*c)", y, y, scs, 0.0, 0.0, -4*pi, +4*pi);
}

cmd apodizace_CB(int slider_id, double apodiz)
{

double souradnice_x1[2*NBR_PTS + 1];
double souradnice_y1[2*NBR_PTS + 1];
double souradnice_x2[2*NBR_PTS + 1];
double souradnice_y2[2*NBR_PTS + 1];



double x; //souradnice v obrazove rovine
double pom1,pom2,pom3;
int i;

/*Vypocet rozdeleni intenzity obrazu bodu pro kruhovou pupilu v definovanych bodech */ 
	for (i = -NBR_PTS; i <= NBR_PTS; i++) //vypocet v 2*NBR_PTS bodech
	{
	  x = 4*pi*i/NBR_PTS; //zakresleni probehne v rozsahu -4*pi az 4*pi
	  souradnice_x1[i + NBR_PTS] = x;  //souradnice x se ulozi do pole
 
	  if(x!=0)
	    souradnice_y1[i + NBR_PTS] = (4*j1(x)*j1(x))/(x*x); //vypoctene hodnoty intenzity se ulozi do pole souradnice_y1[]
	  else
	    souradnice_y1[i+NBR_PTS]=1;
	}
	
/*Vypocet rozdeleni intenzity obrazu bodu pro pupilu jejiz propustnost je promenna dle vztahu g(q)= 1-konst*q    */ 
	for (i = -NBR_PTS; i <= NBR_PTS; i++)
	{
	  x = 4*pi*i/NBR_PTS; //zakresleni probehne v rozsahu -4*pi az 4*pi
	  souradnice_x2[i + NBR_PTS] = x; //souradnice x se ulozi do pole

	
	  pom1=j1(x);
	  pom2=j0(x);

	  if(x!=0)
	    { 
 pom3=(2.0/(x*x))*((x*(1.0-apodiz)+4*apodiz/x)*pom1-2*apodiz*pom2);
	      souradnice_y2[i + NBR_PTS] = pom3*pom3; //vypoctene hodnoty intenzity se ulozi do pole souradnice_y2[]

	    }
	      else 
		souradnice_y2[i+NBR_PTS]=1;
	    
	}
	
	// nasleduje kresleni
	graphwin_open(SLIDERBASEWIN);
	eperase;
	gplot("2", 2*NBR_PTS + 1,
	      souradnice_x1 , souradnice_y1 , "kruhova pupila",
	      souradnice_x2, souradnice_y2," mezikruzi",
	      "relativni intenzita", "x=pi*y'/(lambda*c)", y, y, scs, 0.0, 0.0, -4*pi, +4*pi);
}