/*
 * Calculate field of view for lenses.
 *
 * Greg Lehey, 28 March 2010.
 *
 * $Id: fov.c,v 1.17 2020/12/15 03:49:58 grog Exp $
 */

#include <math.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <sys/file.h>
#include <sys/stat.h>
#include <sys/types.h>
#include <unistd.h>
#include "photos.h"

void usage (char *me)
{
  fprintf (stderr,
           "Usage: %s [-a fov][-d distance] [-f] [-m magnification] [-O overlap] [-s sensor] [focal-length ...]\n\n"
           "sensor can be:\n"
           "     F\tFull frame\n"
           "     H\tAPS H\n"
           "     C\tAPS C (Nikon)\n"
           "     Canon\tAPS C (Canon)\n"
           "     S\tFoveon (Sigma)\n"
	   "-f specifies fisheye\n"
           "fov is the field of view in degrees, overrides focal length\n"
           "distance is measured from the entrance pupil\n"
           "magnification is the magnification for macros\n"
           "overlap is percentage overlap for panoramas (default 30%%)\n",
           me );
  exit (1);
}

int main (int argc, char *argv [])
{
  int i;                                        /* index in argv */
  int fisheye = 0;				/* not a fisheye projection */
  struct sensor *s = &fourthirds;
  float f;                                      /* focal length */
  float overlap;                                /* percentage overlap for panoramas */
  float distance = 0;				/* focus distance.  0 means infinity :-) */
  float magnification = 0;			/* magnificaion.  0 means use distance */
  float dfov = 0.0;				/* field of view (degrees), also flag */
  float fov;					/* field of view in radians */

  float hangle;                                 /* horizontal FoV */
  float vangle;                                 /* vertical FoV */
  float dangle;					/* diagonal FoV */

  i = 1;                                        /* first arg to process */
  overlap = 0.3;                                /* overlap 30% by default */
  fullframe.diagonal = sqrt (fullframe.width * fullframe.width + fullframe.height * fullframe.height);
  fullframe.aspect_ratio = fullframe.width / fullframe.height;

  if (argc < 2)                                 /* no parameters, */
    usage (argv [0]);

  for (i = 1; i < argc; i++)
    {
    if (! strcmp (argv [i], "-a"))		/* -a: by field of view */
      {
      i++;					/* point to arg */
      if (i == argc)
        usage (argv [0]);
      dfov = atof (argv [i]);			/* field of view in degrees */
      fov = dfov / degrees;			/* and in radians */
      }
    else if (! strcmp (argv [i], "-f"))		/* -d: distance from entrance pupil */
      fisheye = 1;				/* fisheye view */
    else if (! strcmp (argv [i], "-d"))		/* -d: distance from entrance pupil */
      {
      i++;					/* point to arg */
      if (i == argc)
        usage (argv [0]);
      distance = atof (argv [i]);
      }
    else if (! strcmp (argv [i], "-m"))		/* -m: magnification */
      {
      i++;					/* point to arg */
      if (i == argc)
        usage (argv [0]);
      magnification = atof (argv [i]);
      }
    else if (! strcmp (argv [i], "-O"))		/* -O: overlap */
      {
      i++;					/* point to arg */
      if (i == argc)
        usage (argv [0]);
      overlap = atof (argv [i]) * 0.01;
      }
    else if (! strcmp (argv [i], "-s"))
      {
      i++;					/* point to arg */
      if (i == argc)
        usage (argv [0]);
      else if (strcmp (argv [i], "F") == 0)
        s = &fullframe;
      else if (strcmp (argv [i], "H") == 0)
        s = &aps_h;
      else if (strcmp (argv [i], "C") == 0)
        s = &aps_c;
      else if (strcmp (argv [i], "Canon") == 0)
        s = &aps_canon;
      else if (strcmp (argv [i], "S") == 0)
        s = &foveon;
      else
        usage (argv [0]);
      }
    else if (argv [i] [0] == '-')
      usage (argv [0]);
    else
      break;
    }

  s->diagonal = sqrt (s->width * s->width + s->height * s->height);
  s->aspect_ratio = s->width / s->height;
  printf ("Sensor details:\n"
          "Width:\t\t\t\t%5.1f mm\n"
          "Height:\t\t\t\t%5.1f mm\n"
          "Diagonal:\t\t\t%5.1f mm\n"
          "Area:\t\t\t\t%5.1f mm²\n"
          "Width ratio to full frame:\t%7.3f\n"
          "Height ratio to full frame:\t%7.3f\n"
          "Diagonal ratio to full frame:\t%7.3f\n"
          "Area ratio to full frame:\t%7.3f\n\n",
          s->width,
          s->height,
          s->diagonal,
          s->width * s->height,
          fullframe.width / s->width,
          fullframe.height / s->height,
          fullframe.diagonal / s->diagonal,
          fullframe.width * fullframe.height / (s->width * s->height) );

  if (dfov)					/* we want focal lengths for this fov */
    {
    /*
     * At least for the time being, fov can be horizontal, vertical and
     * diagonal.  Do all three.
     */
    float df;					/* diagonal "focal length" */
    float hf;					/* horizontal "focal length" */
    float vf;					/* vertical "focal length" */
    if (fisheye)
      {
      df = s->diagonal / fov;
      hf = s->width / fov;
      vf = s->height / fov;
      }
    else					/* rectilinear */
      {
      float mytan = 2 * tan (fov / 2);

      if (dfov > 179)
	{
	fprintf (stderr, "Field of view too large for rectilinear projection (%f > 179)\n", dfov);
	exit (1);
	}
      df = s->diagonal / mytan;
      hf = s->width / mytan;
      vf = s->height / mytan;
      }

    printf ("Field of view:\t\t\t%6.2f°\n"
	    "Diagonal focal length:\t\t%6.2f mm\n"
	    "Horizontal focal length:\t%6.2f mm\n"
	    "Vertical focal length:\t\t%6.2f mm\n",
	    dfov,
	    df,
	    hf,
	    vf );
    }
  else
    for (; i < argc; i++)
      {
      f = atof (argv [i]);
      if (fisheye)
      /*
       * What projection is our fisheye lens?  I'm assuming equidistant here,
       * but http://bobatkins.com/photography/technical/field_of_view.html comes
       * up with formulae for other projections:
       *
       *    FOV (equisolid fisheye) = 4 * arcsin (frame size/(focal length * 4))
       *    FOV  (orthogonal fisheye) = 2 * arcsin (frame size/(focal length *2)
       *    FOV (stereographic fisheye) = 4 * arctan (frame size/(focal length * 4))
       */
	{
	hangle = s->width / f * degrees;
	vangle = s->height / f * degrees;
	dangle = s->diagonal / f * degrees;
	}
      else
	{
	hangle = atan (s->width / (f * 2.0)) * 2.0 * degrees;
	vangle = atan (s->height / (f * 2.0)) * 2.0 * degrees;
	dangle = atan (s->diagonal / (f * 2.0)) * 2.0 * degrees;
	}

      printf ("Focal length:\t\t\t%6.2f mm\n"
	      "Horizontal FOV:\t\t\t%6.2f°\n"
	      "Diagonal FOV:\t\t\t%6.2f°\n"
	      "Vertical FOV:\t\t\t%6.2f°\n",
	      f,
	      hangle,
	      dangle,
	      vangle );

      if (magnification)				/* override distance */
	distance = f / magnification - f;
      if (distance)
	/*
	 * This is the same as above, except that v (distance from sensor) takes
	 * place of f.
	 */
	{
	float v = 1 / (1 / f - 1 / (distance * 1000)); /* distance from sensor */

	if (fisheye)
	  {
	  hangle = s->width / (v * 2.0) * 2.0 * degrees;
	  vangle = s->height / (v * 2.0) * 2.0 * degrees;
	  dangle = s->diagonal / (v * 2.0) * 2.0 * degrees;
	  }
	else
	  {
	  hangle = atan (s->width / (v * 2.0)) * 2.0 * degrees;
	  vangle = atan (s->height / (v * 2.0)) * 2.0 * degrees;
	  dangle = atan (s->diagonal / (v * 2.0)) * 2.0 * degrees;
	  }

	printf ("At %1.2f m:\n"
		"Effective focal length:\t\t%6.2f mm\n"
		"Horizontal FOV:\t\t\t%6.2f°\n"
		"Diagonal FOV:\t\t\t%6.2f°\n"
		"Vertical FOV:\t\t\t%6.2f°\n",
		distance,
		v,
		hangle,
		dangle,
		vangle );
	}

      printf ("Panorama steps (%1.0f%% overlap):\n"
	      "Horizontal:\t\t\t%3.0f°\n"
	      "Vertical:\t\t\t%3.0f°\n\n",
	      overlap * 100,
	      hangle * (1 - overlap),
	      vangle * (1 - overlap) );
      }

  return 0;
}