#include <windows.h>
#include <math.h>
#include <string.h>
#include "usersurf.h"

/*
Written by Kenneth E. Moore
June 1, 1998
*/

int __declspec(dllexport) APIENTRY UserDefinedSurface(USER_DATA *UD, FIXED_DATA *FD);

/* a generic Snells law refraction routine */
int Refract(double thisn, double nextn, double *l, double *m, double *n, double ln, double mn, double nn);

BOOL WINAPI DllMain (HANDLE hInst, ULONG ul_reason_for_call, LPVOID lpReserved)
	{
   return TRUE;
   }

/*

This DLL models a plane surface, with an additional wavelength dependent propagation
hidden within the surface to model systems with independent focus in each wavelength.

*/

int __declspec(dllexport) APIENTRY UserDefinedSurface(USER_DATA *UD, FIXED_DATA *FD)
	{
   int i;
   double t;
   switch(FD->type)
   	{
      case 0:
      	/* ZEMAX is requesting general information about the surface */
         switch(FD->numb)
         	{
            case 0:
            	/* ZEMAX wants to know the name of the surface */
		         /* do not exceed 12 characters */
		         strcpy(UD->string,"Offset");
               break;
            case 1:
            	/* ZEMAX wants to know if this surface is rotationally symmetric */
               /* it is, so return any character in the string; otherwise, return a null string */
            	strcpy(UD->string, "1");
               break;
            case 2:
            	/* ZEMAX wants to know if this surface is a gradient index media */
               /* it is not, so return a null string */
            	UD->string[0] = '\0';
            	break;
            }
         break;
      case 1:
      	/* ZEMAX is requesting the names of the parameter columns */
         /* the value FD->numb will indicate which value ZEMAX wants. */
         /* they are all "Unused" for this surface type */
         /* returning a null string indicates that the parameter is unused. */
         switch(FD->numb)
         	{
            default:
            	UD->string[0] = '\0';
            	break;
            }
      	break;
      case 2:
      	/* ZEMAX is requesting the names of the extra data columns */
         /* the value FD->numb will indicate which value ZEMAX wants. */
         /* they are all "Unused" for this surface type */
         /* returning a null string indicates that the extradata value is unused. */
         switch(FD->numb)
         	{
            case 1:
		         strcpy(UD->string,"Wv 1 Offset");
               break;
            case 2:
		         strcpy(UD->string,"Wv 2 Offset");
               break;
            case 3:
		         strcpy(UD->string,"Wv 3 Offset");
               break;
            case 4:
		         strcpy(UD->string,"Wv 4 Offset");
               break;
            case 5:
		         strcpy(UD->string,"Wv 5 Offset");
               break;
            case 6:
		         strcpy(UD->string,"Wv 6 Offset");
               break;
            case 7:
		         strcpy(UD->string,"Wv 7 Offset");
               break;
            case 8:
		         strcpy(UD->string,"Wv 8 Offset");
               break;
            case 9:
		         strcpy(UD->string,"Wv 9 Offset");
               break;
            case 10:
		         strcpy(UD->string,"Wv 10 Offset");
               break;
            case 11:
		         strcpy(UD->string,"Wv 11 Offset");
               break;
            case 12:
		         strcpy(UD->string,"Wv 12 Offset");
               break;
            default:
            	UD->string[0] = '\0';
            	break;
            }
      	break;
      case 3:
      	/* ZEMAX wants to know the sag of the surface */
         /* it is always zero */
         UD->sag1 = 0.0;
         UD->sag2 = 0.0;
      	break;
      case 4:
      	/* ZEMAX wants a paraxial ray trace to this surface */
         /* x, y, z, and the optical path are unaffected, at least for this surface type */
         /* for paraxial ray tracing, the return z coordinate should always be zero. */
         /* paraxial surfaces are always planes with the following normals */

         UD->ln =  0.0;
         UD->mn =  0.0;
         UD->nn = -1.0;

         if ((UD->n) != 0.0)
         	{
            (UD->l) = (UD->l)/(UD->n);
            (UD->m) = (UD->m)/(UD->n);

            (UD->l) = (FD->n1*(UD->l))/(FD->n2);
            (UD->m) = (FD->n1*(UD->m))/(FD->n2);

            /* normalize */
            (UD->n) = sqrt(1/(1 + (UD->l)*(UD->l) + (UD->m)*(UD->m) ) );
            /* de-paraxialize */
            (UD->l) = (UD->l)*(UD->n);
            (UD->m) = (UD->m)*(UD->n);

   	      /* now propagate the offset distance */
      	   t = FD->xdata[FD->wave] / (UD->n);
				(UD->x) = (UD->l) * t + (UD->x);
				(UD->y) = (UD->m) * t + (UD->y);
            }

         break;
      case 5:
      	/* ZEMAX wants a real ray trace to this surface */

         if ((UD->n) != 0.0)
         	{
	         /* propagate the offset distance */
   	      t = FD->xdata[FD->wave] / (UD->n);
				(UD->x) = (UD->l) * t + (UD->x);
				(UD->y) = (UD->m) * t + (UD->y);
            }
			UD->path = t;

         UD->ln =  0.0;
         UD->mn =  0.0;
         UD->nn = -1.0;
         if (Refract(FD->n1, FD->n2, &UD->l, &UD->m, &UD->n, UD->ln, UD->mn, UD->nn)) return(-FD->surf);
         break;
      case 6:
      	/* ZEMAX wants the index, dn/dx, dn/dy, and dn/dz at the given x, y, z. */
         /* This is only required for gradient index surfaces, so return dummy values */
         UD->index = FD->n2;
         UD->dndx = 0.0;
         UD->dndy = 0.0;
         UD->dndz = 0.0;
      	break;
      case 7:
      	/* ZEMAX wants the "safe" data. */
         /* this is used by ZEMAX to set the initial values for all parameters and extra data */
         /* when the user first changes to this surface type. */
         /* this is the only time the DLL should modify the data in the FIXED_DATA FD structure */
         for (i = 1; i <= 8; i++) FD->param[i] = 0.0;
         for (i = 1; i <= 200; i++) FD->xdata[i] = 0.0;
         break;
      }
   return 0;
   }

int Refract(double thisn, double nextn, double *l, double *m, double *n, double ln, double mn, double nn)
{
double nr, cosi, cosi2, rad, cosr, gamma;
if (thisn != nextn)
	{
	nr = thisn / nextn;
	cosi = fabs((*l) * ln + (*m) * mn + (*n) * nn);
	cosi2 = cosi * cosi;
	if (cosi2 > 1) cosi2 = 1;
	rad = 1 - ((1 - cosi2) * (nr * nr));
	if (rad < 0) return(-1);
	cosr = sqrt(rad);
	gamma = nr * cosi - cosr;
	(*l) = (nr * (*l)) + (gamma * ln);
	(*m) = (nr * (*m)) + (gamma * mn);
	(*n) = (nr * (*n)) + (gamma * nn);
	}
return 0;
}

