--
-- To execute these examples you should login as user "root" with:
--
-- mysql --local-infile=1 -u root -p
--
-- and type the password when requested.
--
-- In the following examples we will assume that DIF has been properly installed.
-- Remember that comment lines start with -- or # !
-- Also note that functions can only return 1 value!
--
-- Use the "DIF" database, but can use any, e.g. mpe2018db
USE DIF;

--
-- Set geographic coordinate and time zone offset (Garching)
--
-- 48.2618996 , 11.6717692
--
SET @Long := 11.7;
SET @Lat := 483;
SET @TZ_offset := -1;

--
--  M1 (Crab Nebula) coordinates:  RA  = 5.575   (Hours)
--                                 DEC = 22.0167 (Degrees);
--
-- To convert from Hours to Degrees simply multiply by 15.0
--
SET @Ra := 5.575 * 15.0;
SET @Decl := 22.0167;

--
-- This variable will be used to perform circular selection
--
-- Radius in degrees
--
SET @Rad := 10;

--
-- From now on all angles will be in Degrees
--

-- Create the table that will contain the Messier catalog
DROP TABLE IF EXISTS Messier;
CREATE TABLE Messier (
	M int NOT NULL, 
	Type CHAR(2) DEFAULT '**',
	Const CHAR(3) DEFAULT '***', 
	Mag FLOAT, 
	Ra  FLOAT, 
	Decl FLOAT, 
	Dist CHAR(20), 
	App_size CHAR(20) DEFAULT 'unknown', 
        PRIMARY KEY(M),
        INDEX(Type)
);

DESCRIBE Messier;
SELECT * FROM Messier;


-- Load data into the table
LOAD DATA LOCAL INFILE './messier.tsv'  
  INTO TABLE Messier;

SELECT * FROM Messier;   -- note the Ra column values


-- Convert Ra to degrees.
UPDATE Messier SET Ra = Ra * 15.0;
SELECT * FROM Messier;


-- Create the table that will contain the description of each object type
DROP TABLE IF EXISTS TypeDescr;
CREATE TABLE TypeDescr(
      Type CHAR(2), 
      Descr CHAR(20), 
      PRIMARY KEY(Type)
);
DESCRIBE TypeDescr;

-- Load data into the table
LOAD DATA LOCAL INFILE './messier_legend.tsv'  
  INTO TABLE TypeDescr;

SELECT * FROM TypeDescr;



-- A simple join to show the long description of each object type
SELECT Messier.*, TypeDescr.Descr 
FROM Messier LEFT JOIN TypeDescr
ON Messier.Type = TypeDescr.Type;



-- From here we assume you have DIF installed!
-- To index the table, from the shell give the command:
-- To index the table, from the shell give the command:
--
--   dif --index-htm DIF Messier 6 Ra Decl
--
SELECT *
FROM Messier_htm_6
WHERE DIF_Circle(@Ra, @Decl, @Rad);



-- Examples of useful functions and procedures

delimiter //
--
-- Compute distance between two points on a sphere
--
--  Note that DIF already has a function to calculate spherical distance: Sphedist
--  This function is just for demonstration.
--
-- a1/b1: long/lat of first point (DEG)
-- a2/b2: long/lat of second point (DEG)
--
-- Returns:
--    distance in Degrees
--
DROP FUNCTION IF EXISTS skydist//
CREATE FUNCTION skydist(a1 DOUBLE, b1 DOUBLE, a2 DOUBLE, b2 DOUBLE) RETURNS DOUBLE
BEGIN
  DECLARE deg2rad, a1r, b1r, a2r, b2r, radiff, cosdis, dist DOUBLE;

  SET deg2rad = 0.017453292519943;
  SET a1r = a1*deg2rad;
  SET b1r = b1*deg2rad;
  SET a2r = a2*deg2rad;
  SET b2r = b2*deg2rad;
  SET radiff = ABS(a2r-a1r);
  IF (radiff > pi()) THEN
    SET radiff = pi() * 2e0 - radiff;
  END IF;
  SET cosdis = sin(b1r)*sin(b2r) + cos(b1r)*cos(b2r)*cos(radiff);
  SET dist = acos(cosdis)/deg2rad;

  RETURN dist;
END//


-- From year, month, day, hour to Julian day
DROP FUNCTION IF EXISTS jdcnv//
CREATE FUNCTION jdcnv(yr INT, mn INT, day INT, hr DOUBLE) RETURNS DOUBLE
BEGIN
  DECLARE L DOUBLE;
  DECLARE j BIGINT;
  DECLARE julian DOUBLE;

  SET L = (mn-14) DIV 12;
  SET j = day - 32075 + 1461*(yr+4800+L) DIV 4 + 
               367*(mn - 2-L*12) DIV 12 - 3*((yr+4900+L) DIV 100) DIV 4;

  SET julian = j + hr/24.0 -0.5;
  RETURN julian;
END//


-- From hour angle and declination to Alt, Az coordinates
DROP PROCEDURE IF EXISTS hadec2altaz//
CREATE PROCEDURE hadec2altaz(IN ha DOUBLE, IN decl DOUBLE, IN lat DOUBLE, 
                             OUT alt DOUBLE, OUT az DOUBLE)
BEGIN
  DECLARE d2r, sh, ch, sd, cd, sl, cl DOUBLE;
  DECLARE x, y, z, r DOUBLE;

  SET d2r = pi()/180.0;
  
  SET sh = sin(ha*d2r);
  SET ch = cos(ha*d2r);
  SET sd = sin(decl*d2r);
  SET cd = cos(decl*d2r);
  SET sl = sin(lat*d2r);
  SET cl = cos(lat*d2r);

  SET x = - ch * cd * sl + sd * cl;
  SET y = - sh * cd;
  SET z = ch * cd * cl + sd * sl;
  SET r = sqrt(x*x + y*y);
  
  SET az = atan(y,x) / d2r;
  SET alt = atan(z,r) / d2r;

  IF (az < 0) THEN
    SET az = az + 360;
  END IF;
END//


-- From hour angle and declination to Altitude
DROP FUNCTION IF EXISTS hadec2alt//
CREATE FUNCTION hadec2alt(ha DOUBLE, decl DOUBLE, lat DOUBLE) RETURNS DOUBLE
DETERMINISTIC
BEGIN
  DECLARE alt, az DOUBLE;
  CALL hadec2altaz(ha, decl, lat, alt, az);

  RETURN alt;
END//


-- From hour angle and declination to Azimuth
DROP FUNCTION IF EXISTS hadec2az//
CREATE FUNCTION hadec2az(ha DOUBLE, decl DOUBLE, lat DOUBLE) RETURNS DOUBLE
DETERMINISTIC
BEGIN
  DECLARE alt, az DOUBLE;
  CALL hadec2altaz(ha, decl, lat, alt, az);

  RETURN az;
END//


-- From Julian day to local sidereal time
DROP FUNCTION IF EXISTS jd2st//
CREATE FUNCTION jd2st(jd DOUBLE, lng DOUBLE) RETURNS DOUBLE
DETERMINISTIC
BEGIN
  DECLARE c0, c1, c2, c3 DOUBLE;
  DECLARE jd2000, lst DOUBLE;
  DECLARE t, t0, theta DOUBLE;

  SET c0 = 280.46061837;
  SET c1 = 360.98564736629;
  SET c2 = 0.000387933;
  SET c3 = 38710000.0;
 
  SET jd2000 = 2451545.0;

  SET t0 = jd - jd2000;
  SET t = t0/36525.0;
  
  SET theta = c0 + (c1 * t0) + t*t*(c2 - t / c3 );
  
  SET lst = (theta + lng)/15.0;

  IF (lst < 0) THEN
    SET lst = 24.0 + (lst % 24.0);
  END IF;
  SET lst = lst % 24.0;

  RETURN lst;
END//


-- From local DATETIME and Time zone to Julian day
DROP FUNCTION IF EXISTS Julian//
CREATE FUNCTION Julian(Time DATETIME, TZ_offset INT) RETURNS DOUBLE
DETERMINISTIC
BEGIN
  DECLARE Yr, Mn, Dy, Hr, Mi, Se INT;
  DECLARE Hr_fr, jd DOUBLE;

  SET Yr = year(Time);
  SET Mn = month(Time);
  SET Dy = day(Time);
  SET Hr = hour(Time);
  SET Mi = minute(Time);
  SET Se = second(Time);
  SET Hr_fr =(Hr+TZ_offset) + Mi/60.0 + Se/3600.0;
  
  SET jd = jdcnv(Yr, Mn, Dy, Hr_fr);
  RETURN jd;
END//


-- From local DATETIME and Time zone to local sidereal time
DROP FUNCTION IF EXISTS Sidereal//
CREATE FUNCTION Sidereal(Time DATETIME, TZ_offset INT, Lng DOUBLE) RETURNS DOUBLE
DETERMINISTIC
BEGIN
  DECLARE jd, st DOUBLE;

  SET jd = Julian(Time, TZ_offset);
  SET st = jd2st(jd, Lng);
  RETURN st;
END//


-- From equatorial Ra, Dec to Altitude
DROP FUNCTION IF EXISTS radec2alt//
CREATE FUNCTION radec2alt(Ra DOUBLE, Decl DOUBLE, Lat DOUBLE, Lng DOUBLE, 
                          TZ_offset INT, Time DATETIME)
RETURNS DOUBLE
BEGIN
  DECLARE st, ha, alt DOUBLE;

  SET st = Sidereal(Time, TZ_offset, Lng) * 15.0;
  SET ha = st - Ra;
  SET alt = hadec2alt(ha, Decl, Lat);
  RETURN alt;
END//


-- From equatorial Ra, Dec to Azimuth
DROP FUNCTION IF EXISTS radec2az//
CREATE FUNCTION radec2az(Ra DOUBLE, Decl DOUBLE, Lat DOUBLE, Lng DOUBLE, 
                         TZ_offset INT, Time DATETIME)
RETURNS DOUBLE
BEGIN
  DECLARE st, ha, az DOUBLE;

  SET st = Sidereal(Time, TZ_offset, Lng) * 15.0;
  SET ha = st - Ra;
  SET az = hadec2az(ha, Decl, Lat);
  RETURN az;
END//
delimiter ;



--
-- Example of usage of jdcnv and jd2st
--
SELECT  @Time:= now() AS Time, 
	@Yr:=year(@Time) AS Year, 
	@Mn:=month(@Time) AS Month, 
	@Dy:=day(@Time) AS Day,
	@Hr:=hour(@Time) AS Hour,
	@Mi:=minute(@Time) AS Minute,
	@Se:=second(@Time) AS Second,
	@Hr_fr:=(@Hr+@TZ_offset) + @Mi/60.0 + @Se/3600.0 AS Hr_fr,
	@Jd:=jdcnv(@Yr, @Mn, @Dy, @Hr_fr) AS Julian,
	@Si:=jd2st(@Jd, @Long) AS Sidereal;

-- Do the same with
SELECT Julian(now(), @TZ_offset), Sidereal(now(), @TZ_offset, @Long);

--
-- Example of usage of hadec2alt, hadec2az
-- (coordinate conversion from Ra/Dec to Alt/Az)
--
-- Compute sidereal time and convert to Degrees
--
SET @Si := Sidereal(now(), @TZ_offset, @Long) * 15.0;

--
-- Compute hour angle
--
SET @Hour_angle := @Si - @Ra;

SELECT @Ra, @Decl, @Hour_angle, 
       hadec2alt(@Hour_angle, @Decl, @Lat) AS Alt, 
       hadec2az( @Hour_angle, @Decl, @Lat) AS Az;

-- Do the same with
SELECT @Ra, @Decl, @Hour_angle,
       radec2alt(@Ra, @Decl, @Lat, @Long, @TZ_offset, now()) AS Alt,
       radec2az( @Ra, @Decl, @Lat, @Long, @TZ_offset, now()) AS Az;


--
-- Note: variables (like @Lat, @Long, @TZ_offset) cannot be used inside
-- views.
--
DROP VIEW IF EXISTS Messier_view;
CREATE VIEW Messier_view AS
  SELECT Messier.M, 
         TypeDescr.Descr AS Descr,
         Messier.Const,
         Messier.Mag,
         Messier.Ra,
         Messier.Decl,
         radec2alt(Messier.Ra, Messier.Decl, 38, 13, -1, now()) AS Alt,
         radec2az( Messier.Ra, Messier.Decl, 38, 13, -1, now()) AS Az,
         Messier.Dist,
         Messier.App_size
  FROM Messier LEFT JOIN TypeDescr
  ON Messier.Type = TypeDescr.Type
  ORDER BY Descr, M;

DESCRIBE Messier_view;
SELECT * FROM Messier_view;


-- *******************************************************************

-- This part deals with the "manual" management of DIF to index tables.
-- You should use the script "dif" to do this, but if you want to practice
-- with DB management you can try these examples.
--
-- Set up the DIF_Rect, DIF_Circle, and HTMLookup functions.
-- 
-- NOTE 1: this require that DIF has been installed and that it is placed
-- in the 'plugin_dir' of MySQL.
--

show variables like 'plugin%';
+---------------+------------------------------+
| Variable_name | Value                        |
+---------------+------------------------------+
| plugin_dir    | /usr/local/mysql/lib/plugin/ |
+---------------+------------------------------+

-- Then from the shell: ls /usr/local/mysql/lib/plugin/udf_dif.so
-- 
-- NOTE 2: after DIF has been installed issue 'ldconfig' from a root
-- shell if you have troubles seeing the functions.
-- 

--
-- DIF_Circle(Ra_deg DOUBLE, Dec_deg DOUBLE, Rad_arcmin DOUBLE)
--
-- Generate htmIDs for those pixels in a circle centered at Ra_deg, Dec_deg
-- (deg) with a radius of Rad_arcmin (arcmin).
--
-- Returns:
--  how many pixels fall into the circle
--
DROP FUNCTION IF EXISTS DIF_Circle;
CREATE FUNCTION DIF_Circle RETURNS INTEGER SONAME 'ha_dif.so';

--
-- HTMLookup(depth INT, Ra_deg DOUBLE, Dec_deg DOUBLE)

-- Lookup the htmID, at a given depth, of the pixel which contains the point
-- with with coordinates Ra_deg and Dec_deg.
--
-- Returns:
--  the HTM id of the pixel at the requested depth
--
DROP FUNCTION IF EXISTS HTMLookup;
CREATE FUNCTION HTMLookup RETURNS INTEGER SONAME 'ha_dif.so';


-- This is to query a rectangular region - see DIF manual.
DROP FUNCTION IF EXISTS DIF_Rect;
CREATE FUNCTION DIF_Rect RETURNS INTEGER SONAME 'ha_dif.so';



--
-- Add a column and an index for the HTM id (depth 6) to the Messier table
--
ALTER TABLE Messier ADD COLUMN htmID_6 SMALLINT UNSIGNED, ADD INDEX(htmID_6);

DESCRIBE Messier;

-- Populate the new column with htmIDs
UPDATE Messier SET htmID_6 = HTMLookup(6, Ra, Decl);
SELECT * FROM Messier;


--
-- Set up a trigger to automatically update the htmID_6 field when a new
-- record is inserted
-- Note that dif might have already has installed triggers for INSERT and UPDATE,
-- so you need to remove them before you can create a new one
--

show triggers;

-- eventually
drop trigger trig_ins_Messier;

delimiter //
DROP TRIGGER trig_ins_Messier//
CREATE TRIGGER trig_ins_Messier BEFORE INSERT ON Messier FOR EACH ROW
BEGIN
  SET NEW.htmID_6 = HTMLookup(6, NEW.Ra, NEW.Decl);
END//

--
-- Do the same for updated records
--
DROP TRIGGER trig_upd_Messier//
CREATE TRIGGER trig_upd_Messier BEFORE UPDATE ON Messier FOR EACH ROW
BEGIN
  SET NEW.htmID_6 = HTMLookup(6, NEW.Ra, NEW.Decl);
END//
delimiter ;
