using System;
using System.IO;
using System.Drawing;
using System.Drawing.Imaging;
using System.Drawing.Drawing2D;
using System.Data;
using System.Data.SqlClient;
using System.Data.SqlTypes;

namespace Geom
{

	/// <summary>
	/// Summary description for FrameWCS.
	/// </summary>
	/// <summary>
	/// FrameWCS carries the transformation of ra/dec to xy and back.
	/// Created by Alex Szalay, 2004-08-16, from the SDSS specific Frame class. 
	/// This class is working on WCS entirely. 
	/// </summary>
	public class FrameWCS
	{
		public static string Revision
		{
			get 
			{
				return "$Revision: 1.1.1.1 $";
			}
		}
		//--------------------------------------
		// these values are specific to NDWFS	
		//--------------------------------------
        public static float     FrameWidth; 			            // pixel
		public static float		FrameHeight;			            // pixel
		public static float		PixelsPerDegree; 		            // pixel
		public static float		ClipX			= 0.0F;				// clip from the left edge
		public static float		ClipY			= 0.0F;				// clip from the top edge
		public static float		PixX			= 0.5F;				// WCS reference within pixel
		public static float		PixY			= 0.5F;
		public static string	SurveyName		= "";   			// survey name
		//--------------------------------------------------------------------------------------------------
		// these are generic, derived from the ones above the line
		//--------------------------------------------------------------------------------------------------
		public static double	PixelsPerArcMin;		// derived from above
		public static double	FrameHalfDiag;	
		private static double	D2R = Math.PI/180.0;							// degrees to radians		
		//--------------------------------------------------------------------------------------------------
		public  double ra,dec;							// point (ra, dec) in J2000 degrees
		public  double scale;							// image scale 
		public  string info;							// info about the field for debugging
		public  string fieldid;							// string representation of the fieldid
		public  Matrix m;
		public  Image  tile;

		private double crpix1, crpix2, crval1, crval2;	// the WCS
		private double cd1_1, cd1_2, cd2_1, cd2_2;		// the WCS
		private double qd1_1, qd1_2, qd2_1, qd2_2;		// inverse transformation matrix
		private float  x,y;								// pixel offsets from the fiducial point
		private double[] n,w,u;							// the three normal vectors
        private string gridMask;                        // GridMask for Galex cutout.

		//-------------------------------------------------------------------------------

		// -----------------------------------------------------------------------------
		private static byte		_sflag		= 1;									//  specObject
		private static byte		_pflag		= 2;									//  photoObject
		private static byte		_tflag		= 4;									//  targetObject

		// ----------------------------------------------------------------------------

        // -----------------------------------------------------------------------------
        private PointF cOffset;							// reference point offset
        // -----------------------------------------------------------------------------

        public static void setFrameWidth (float frameWidth)
        {
            FrameWCS.FrameWidth = frameWidth;
        }
        public static void setFrameHeight(float frameHeight)
        {
            FrameWCS.FrameHeight = frameHeight;
        }
        public static void setPixelsPerDegree(float ppd)
        {
            FrameWCS.PixelsPerDegree = ppd;
        }
        public static void setPixlesPerArcMin(float ppam)
        {
            FrameWCS.PixelsPerArcMin = ppam;
        }

        public static double GetFrameHalfDiag ()
        {
            FrameHalfDiag = Math.Sqrt(FrameWidth*FrameWidth+FrameHeight*FrameHeight) / PixelsPerArcMin / 2.0;
            return FrameWCS.FrameHalfDiag;
        }

		/// <summary>
		/// GetFrameSQL returns the SQL string required to fetch Frame
		/// </summary>
		/// <param name="raCenter">right ascension of center point</param>
		/// <param name="decCenter">declination of center point</param>
		/// <param name="fradius">frame radius in arcminutes</param>
		/// <param name="zoom">integer zoom factor (0 is unscaled)</param>
		public static String GetFrameSQL(double raCenter, double decCenter, 
			double fradius, int zoom, bool drawFrames )  
		{
			int zoom10x = (zoom>3)? 30: zoom * 10;		// for DB lookup, clipped at 30
			String sQ = "";

			sQ +=	 "select n.crpix1, n.crpix2, n.crval1, n.crval2,\n"
				+"  n.cd1_1, n.cd1_2, n.cd2_1, n.cd2_2" 
				+", n.fieldID, cast(n.fieldID as varchar(8)) as info";

			//------------------------------------------------
			// only read the JPG bytes if we draw the frames
			//------------------------------------------------
			if (drawFrames) sQ += ", f.img as tile, f.GridMask";

            sQ += "\n from  dbo.fGetNearbyGalexFrameEq ("
				+raCenter.ToString("F6")+","+decCenter.ToString("F6")+","
				+fradius.ToString("F6")+","+zoom10x.ToString()+",'z') n"
				+", FrameWCS_test f \n"
				+"where f.fieldID = n.fieldID and f.zoom ="
				+zoom10x.ToString()+"\n";
				
			return (String) sQ;
		}



		//----------------------------------
		// the rest is not survey specific
		//----------------------------------
		
		public FrameWCS()
        { cOffset = new PointF(0, 0); }

		public FrameWCS(SqlDataReader reader, float sc, bool drawFrames)
		{
			crpix1  = Convert.ToDouble(reader[0]);	// crpix1
			crpix2  = Convert.ToDouble(reader[1]);	// crpix2
			crval1  = Convert.ToDouble(reader[2]);	// crval1
			crval2  = Convert.ToDouble(reader[3]);	// crval2
			cd1_1   = Convert.ToDouble(reader[4]);	// cd1_1
			cd1_2   = Convert.ToDouble(reader[5]);	// cd1_2
			cd2_1   = Convert.ToDouble(reader[6]);	// cd2_1
			cd2_2   = Convert.ToDouble(reader[7]);	// cd2_2
			fieldid = Convert.ToString(reader[8]);	// fieldid
			info    = Convert.ToString(reader[9]);	// info
            

			scale   = sc;							// scale

            #region IF reading img from the Database
            ///<comment>
            /// This is the code if you are reading
            /// the image from the database.
            ///</comment>
            #endregion

            //tile = null;
            //if (drawFrames)
            //    tile = Image.FromStream(new MemoryStream((Byte[])reader[10]));

            Image tmpImage = null;
            tile = null;
            if (drawFrames)
            {
                tmpImage = Image.FromStream(new MemoryStream((Byte[])reader[10]));
                gridMask = Convert.ToString(reader[11]);// GridMask
                tile = ApplyGridMask(tmpImage, gridMask);
            }

            # region IF reading img from URL
            ///<comment>
            /// If you want to read the img from a URL, uncomment everything between the region
            /// and comment out everything that is in the coment above.
            ///</comment>
            //string imgURL = Convert.ToString(reader[10]);
            ////scale = sc;						// scale
            //scale = 1;
            //tile = null;
            //Uri myURL = new Uri(imgURL);
            //tile = Image.FromStream(System.Net.WebRequest.Create(myURL).GetResponse().GetResponseStream());
            #endregion

            // compute the inverse transformation matrix
			double det = (cd1_1*cd2_2-cd1_2*cd2_1);
			qd1_1 =  cd2_2/det;
			qd1_2 = -cd2_1/det;
			qd2_1 = -cd1_2/det;
			qd2_2 =  cd1_1/det;

			// compute the normal vectors
			n   = V3.Normal(crval1,crval2);
			w   = V3.West(crval1,crval2);
			u   = V3.North(crval1,crval2);

            cOffset = new PointF(0, 0);
		}

        public Image ApplyGridMask(Image tmpImage, string gridMask)
        {
            string[] splitArray = gridMask.Split(' ');
            PointF[] mask = new PointF[(int)splitArray.Length/2];
            for (int i = 0; i < splitArray.Length; i=i+2)
            {
                int j = (int)(i / 2);
                mask[j].X = Single.Parse(splitArray[i]);
                //mask[j].Y = tmpImage.Height - Single.Parse(splitArray[i + 1]);
                mask[j].Y = Single.Parse(splitArray[i + 1]);
            }
            tmpImage.Save(@"C:\temp\tmpimage_before_bmp.jpg",ImageFormat.Jpeg);
            Bitmap bmp = new Bitmap(tmpImage.Width, tmpImage.Height, PixelFormat.Format32bppArgb);
            //Bitmap bmp = new Bitmap(tmpImage.Width, tmpImage.Height, PixelFormat.Format24bppRgb);
            bmp.SetResolution(tmpImage.HorizontalResolution, tmpImage.HorizontalResolution);
            GraphicsPath path = new GraphicsPath();
            path.AddLines(mask);
            Graphics grfx = Graphics.FromImage(bmp);
            grfx.SetClip(path);
            grfx.DrawImage(tmpImage, 0, 0);
            bmp.Save(@"C:\temp\tmpimage_after_bmp.jpg", ImageFormat.Jpeg);
            return bmp;
        }

        /// <summary>
        /// For Galex
        /// </summary>
        /// <param name="reader"></param>
        /// <param name="sc"></param>
        public FrameWCS(SqlDataReader reader, float sc, double ra, double dec)
        {
            crpix1 = Convert.ToDouble(reader[0]);	// crpix1
            crpix2 = Convert.ToDouble(reader[1]);	// crpix2
            crval1 = Convert.ToDouble(reader[2]);	// crval1
            crval2 = Convert.ToDouble(reader[3]);	// crval2
            cd1_1 = Convert.ToDouble(reader[4]);	// cd1_1
            cd1_2 = 0.0;	// cd1_2
            cd2_1 = 0.0;	// cd2_1
            cd2_2 = Convert.ToDouble(reader[5]);	// cd2_2
            info = "No Info";
            fieldid = Convert.ToString(reader[6]);	// fieldid
            string imgURL = Convert.ToString(reader[7]);

            //scale = sc;							// scale
            scale = 1;
            tile = null;
            Uri myURL = new Uri(imgURL);
            tile = Image.FromStream(System.Net.WebRequest.Create(myURL).GetResponse().GetResponseStream());

            // compute the inverse transformation matrix
            double det = (cd1_1 * cd2_2 - cd1_2 * cd2_1);
            qd1_1 = cd2_2 / det;
            qd1_2 = -cd2_1 / det;
            qd2_1 = -cd1_2 / det;
            qd2_2 = cd1_1 / det;

            // compute the normal vectors
            n = V3.Normal(crval1, crval2);
            w = V3.West(crval1, crval2);
            u = V3.North(crval1, crval2);

            // Added to test
            cOffset = new PointF(0.0F, 0.0F);
            cOffset = EqToFrame(ra, dec);
        }

		/// <summary>
		/// The copy method copies the values of the coord parameter into the current coordinate object
		/// </summary>
		/// <param name=coord>			
		public void Copy( FrameWCS c) 
		{
			this.crpix1		= c.crpix1;
			this.crpix2		= c.crpix2;
			this.crval1		= c.crval1;
			this.crval2		= c.crval2;
			this.cd1_1		= c.cd1_1;
			this.cd1_2		= c.cd1_2;
			this.cd2_1		= c.cd2_1;
			this.cd2_2		= c.cd2_2;
			this.scale		= c.scale;
			this.info		= c.info;
			this.m			= c.m;
			this.tile		= c.tile;
			
			this.qd1_1		= c.qd1_1;
			this.qd1_2		= c.qd1_2;
			this.qd2_1		= c.qd2_1;
			this.qd2_2		= c.qd2_2;

			this.n			= V3.Copy(c.n);
			this.w			= V3.Copy(c.w);
			this.u			= V3.Copy(c.u);
		}

		/// <summary>
		/// EqtoFrame goes from equatorial coordinates to xy
		/// Sets ra, dec according to the astrometric transformation.
		/// </summary>
		/// <param name="ra">right ascension degrees in J2000</param>
		/// <param name="dec">declination degrees in J2000</param>
		public PointF EqToFrame(double pRa, double pDec) 
		{   
			double xt,yt, xx,yy;
			//-----------------
			// save the angles
			//-----------------
			ra	= pRa;
			dec	= pDec;
			double[] r   = V3.Normal(pRa, pDec);
			//-------------------------------------
			// do the dot products with the vectors
			//-------------------------------------
			double zt  = V3.Dot(r,n)*D2R;
			if (zt==0.0) zt = 0.00000000001;
			xt  = V3.Dot(r,w)/zt;
			yt  = V3.Dot(r,u)/zt;
			
			//----------------------------------------------
			// apply the inverse transform into pixel space
			//----------------------------------------------
			xx = crpix1 + qd1_1*xt + qd1_2*yt;
			yy = crpix2 + qd2_1*xt + qd2_2*yt;

			//---------------------------
			// return pixel coordinates
			//---------------------------
			x = (float)(xx/scale);
			y = (float)(yy/scale);
			PointF p = new PointF(x,y);
			return p;
		} 

		/// <summary>
		/// FrameToEq goes from Frame x,y coordinates and sets ra and dec
		/// </summary>
		/// <param name="x">pixels across (scaled for zoom level)</param>
		/// <param name="y">pixels down (scaled for zoom level)</param>
		public PointEq FrameToEq(float x, float y)  
		{
			//------------------------------------
			// implements the TAN projection only
			//------------------------------------
			// intermediate tangent plane coordinates
			// expressed in radians
			//----------------------------------------
          	double xt  = D2R*(cd1_1*(x-crpix1)+cd1_2*(y-crpix2));
			double yt  = D2R*(cd2_1*(x-crpix1)+cd2_2*(y-crpix2));

			double[] r = V3.Copy(n);			
			for (int i=0;i<3;i++) {
				r[i] += xt*w[i]+yt*u[i];
			}
			PointEq p = V3.ToEq(V3.Normalize(r));
			ra	= p.ra;
			dec = p.dec;
			return p;
		}

		/// <summary>
		/// EqtoFrame goes from equatorial coordinates to xy
		/// Sets ra, dec according to the astrometric transformation.
		/// </summary>
		/// <param name="e">Position as PointEq</param>
		public PointF EqToFrame(PointEq e) 
		{   
			return EqToFrame(e.ra, e.dec);
		}

		/// <summary>
		/// FrameToEq goes from Frame x,y coordinates and sets ra and dec
		/// </summary>
		/// <param name="p">The x,y coordinate as PointF</param>
		public PointEq FrameToEq(PointF p)  
		{
			return FrameToEq(p.X, p.Y);
		}		

		/// <summary>
		/// FieldGeometry. Returns the four corners of the field in 
		/// unscaled pixel coordinates with optional clipping
		///  </summary>
		public static PointF[]	FieldGeometry(bool clip)
		{
			//-------------------------------------------------------------------
			// compute the pixel coordinates of corners with optional clipping
			// add 0.5 to compensate for the IRAF convention
			//-------------------------------------------------------------------
			int clipit	= (clip?1:0);
			float fX1	= (float) (ClipX * clipit + PixX);
			float fY1	= (float) (ClipY * clipit + PixY);
			float fX2	= (float) (FrameWidth  - ClipX * clipit + PixX);
			float fY2	= (float) (FrameHeight - ClipY * clipit + PixY);

			//------------------------------------
			// create array for the four corners
			//------------------------------------
			PointF[] p = new PointF[4];

			p[0] = new PointF(fX1, fY1);
			p[1] = new PointF(fX2, fY1);
			p[2] = new PointF(fX2, fY2);
			p[3] = new PointF(fX1, fY2);

			return p;
		}

		public static byte sflag 
		{
			get 
			{
				return _sflag;
			}
		}

		public static byte pflag 
		{
			get 
			{
				return _pflag;
			}
		}

		public static byte tflag 
		{
			get
			{
				return _tflag ;
			}
		}

	} // end of Frame Class
}


///<summary>
/// Code that was used when the images were not in the database.
/// and we were trying to see the difference between using
/// images of size 3840 X 3840 Vs 512 X 512
///</summary>
//----------------------------------------------------
//        /// <summary>
//        /// For Galex Frames only.
//        /// </summary>
//        /// <param name="raCenter"></param>
//        /// <param name="decCenter"></param>
//        /// <param name="fradius"></param>
//        /// <param name="zoom"></param>
//        /// <returns></returns>
//        public static String GetGalexFrameSQL(double raCenter, double decCenter,
//            double fradius, int zoom, bool drawFrames)
//        {
//            int zoom10x = (zoom>3)? 30: zoom * 10;		// for DB lookup, clipped at 30
//            String sQ = "";

//            sQ += "select n.crpix1, n.crpix2, n.crval1, n.crval2,\n"
//                            + "  n.cd1_1, n.cd1_2, n.cd2_1, n.cd2_2"
//                            + ", n.fieldID, cast(n.fieldID as varchar(8)) as info";

//            //------------------------------------------------
//            // only read the JPG bytes if we draw the frames
//            //------------------------------------------------
//            if (drawFrames) sQ += ", f.img ";

//#region Try with Galex frames (zoom 0) of 3840 X 3840

//            //sQ += "\n from  dbo.fGetNearbyFrameEq ("
//            //    + raCenter.ToString("F6") + "," + decCenter.ToString("F6") + ","
//            //    + fradius.ToString("F6") + "," + zoom10x.ToString() + ",'z') n"
//            //    + ", FrameWCS_test f \n"
//            //    + "where f.fieldID = n.fieldID and f.zoom ="
//            //    + zoom10x.ToString() + "\n";
//#endregion

/////<summary>
///// Trying with the database returning images at zoom 0 of
///// size 480X480 images
//            /// </summary>

//#region Try with galex frames (zoom 0) of 480 X 480

//            sQ += "\n from  dbo.fGetNearbyGalexFrameEq ("
//             + raCenter.ToString("F6") + "," + decCenter.ToString("F6") + ","
//             + fradius.ToString("F6") + "," + zoom10x.ToString() + ",'z') n"
//             + ", galex_cut_images_test_bis f \n"
//             + "where f.fieldID = n.fieldID and f.zoom ="
//             + zoom10x.ToString() + "\n";

//            return (String)sQ;
//#endregion

//            ///<summary>
///// Old stuff for Galex

//            //// Convert the radius in degrees.
//            //fradius = fradius;
//            //double searchRadius = fradius / 60.0;

//            //sQ += "select n.crpix1, n.crpix2, n.crval1, n.crval2,\n"
//            //    + "  n.cdelt1, n.cdelt2, fieldid,"
//            //    + " img";

//            //sQ +="\n from dbo.fGetGalexNearbyFrame (" 
//            //    + raCenter.ToString () + "," + decCenter.ToString () + ","
//            //    + searchRadius.ToString () + "," + zoom10x.ToString () + ")" + " as n" + "\n";

//            //return (String)sQ;
/////</summary>

//        }
//------------------------------------------------------------------------------------------------------