package coordtrans;

import java.util.*;

public class CoordTrans {
    
   /* Transform right acsension and declination in sexadecimal hours and degrees 
    to decimal degrees*/
    public double[] sex2Dec(String ra_hr, String ra_m, String ra_s,
    String de_deg, String de_m, String de_s) throws NumberFormatException
    {
        
        double i_ra_hr = Double.parseDouble(ra_hr);
        double i_ra_m  = Double.parseDouble(ra_m);
        double i_ra_s  = Double.parseDouble(ra_s);
        
        double i_de_deg = Double.parseDouble(de_deg);
        double i_de_m   = Double.parseDouble(de_m);
        double i_de_s   = Double.parseDouble(de_s);
        
        double o_ra_deg;
        double o_de_deg;
        
        o_ra_deg = 15*(i_ra_hr+(i_ra_m+i_ra_s/60)/60);
        o_de_deg = (i_de_deg/Math.abs(i_de_deg))*
                    (Math.abs(i_de_deg)+(i_de_m+i_de_s/60)/60);
        
        double[] output= {o_ra_deg,o_de_deg};
        
        return output;
        
    }
    
   /* Transform right acsension and declination in decimal degrees to 
    sexadecimal hours and degrees*/
    public String[] dec2Sex(String ra, String de) throws NumberFormatException
    {
        
        double i_ra = Double.parseDouble(ra);
        double i_de = Double.parseDouble(de);
      
        int o_ra_hr;
        int o_ra_m;
        double o_ra_s;
        
        int o_de_deg;
        int o_de_m;
        double o_de_s;
        
        o_ra_hr = (int) Math.floor(i_ra/15);
        o_ra_m  = (int) Math.floor(60*(i_ra/15 -o_ra_hr));
        o_ra_s  = 60*((i_ra/15-o_ra_hr)*60 - o_ra_m);
        
        o_de_deg = (int) (Math.abs(i_de)/i_de)*(int) Math.floor(Math.abs(i_de));
        o_de_m   = (int) Math.floor((Math.abs(i_de)-Math.floor(Math.abs(i_de)))*60); 
        o_de_s   = 60*((Math.abs(i_de)-Math.floor(Math.abs(i_de)))*60-o_de_m); 
        
        String o_ra_str = o_ra_hr+" "+o_ra_m+" "+o_ra_s;
        String o_de_str = o_de_deg+" "+o_de_m+" "+o_de_s;
        
        String[] output = {o_ra_str,o_de_str};

        return output;
        
    }    

   /* Precess coordinates from any epoch to another. This calculation 
    is only exact for changing 1950.0 to 2000.0 but should be accurate 
    enough for other years around those dates. Input and output is RA 
    and DEC in decimal degrees.*/
    public double[] precess(double ra, double de, double input_epoch, double output_epoch) 
    {
        double m2 = 1.28076E-2;
        double m1 = 1.28053E-2;
        double n2 = 5.5667E-3;
        double n1 = 5.5674E-3;
        double dy = output_epoch-input_epoch;
        
        double ra_rad = Math.toRadians(ra);
        double de_rad = Math.toRadians(de);
        
        double c = Math.toRadians((m1 + n1*Math.sin(ra_rad)*Math.tan(de_rad))*dy);
        double b = Math.toRadians(n1*Math.cos(ra_rad)*dy);

        double cp = Math.toRadians((m2 + n2*Math.sin(ra_rad+c)*Math.tan(de_rad+b))*dy);
        double bp = Math.toRadians(n2*Math.cos(ra_rad+c)*dy);
        
        double ra_out = ra + Math.toDegrees( (c+cp)/2 );
        double de_out = de + Math.toDegrees( (b+bp)/2 );
        
        double[] output = {ra_out,de_out};
        return output;
        
    }    
    
   /* Converts from Galactic lattitude and longitude (in degrees) to 
   J2000 equatorial coordinates in decimal degrees.*/
    public double[] gal2Equ(double l, double b) 
    {
        double l0    = 32.93;
        double ra0   = 282.86;
        double obliq = 62.87170366;
        double cosobliq = Math.cos(Math.toRadians(obliq));
        double sinobliq = Math.sin(Math.toRadians(obliq));
        double l1    = Math.toRadians(l - l0);
        
        double b_rad = Math.toRadians(b);
        
        double dec_rad = Math.asin(Math.sin(b_rad)*cosobliq + Math.cos(b_rad)*Math.sin(l1)*sinobliq);
        double l2      = Math.acos(Math.cos(b_rad)*Math.cos(l1)/Math.cos(dec_rad));
        
        double ra = 0;
        double sinl2 = (Math.cos(b_rad)*Math.sin(l1)*cosobliq - Math.sin(b_rad)*sinobliq)/Math.cos(dec_rad);
        if(sinl2>0){
            ra = Math.toDegrees(l2)+ra0;
        }else{
            ra = Math.toDegrees(2*Math.PI-l2)+ra0;
        }
        
        if(ra>360){
            ra = ra-360;
        }
        
        double dec = Math.toDegrees(dec_rad);
        
        double[] output = {ra,dec};
        return output;
        
    }    
        

   /* Converts from J2000 equatorial coordinates in decimal degrees
    to Galactic lattitude and longitude (in degrees).*/
    public double[] equ2Gal(double ra, double dec) 
    {
        double l0    = 32.93;
        double ra0   = 282.86;
        double obliq = 62.87170366;
        double cosobliq = Math.cos(Math.toRadians(obliq));
        double sinobliq = Math.sin(Math.toRadians(obliq));
        double l1    = Math.toRadians(ra - ra0);
        
        double dec_rad = Math.toRadians(dec);
        
        double b_rad = Math.asin(Math.sin(dec_rad)*cosobliq - Math.cos(dec_rad)*Math.sin(l1)*sinobliq);
        double l2    = Math.acos(Math.cos(dec_rad)*Math.cos(l1)/Math.cos(b_rad));
        
        double l = 0;
        double sinl2 = (Math.cos(dec_rad)*Math.sin(l1)*cosobliq + Math.sin(dec_rad)*sinobliq)/Math.cos(b_rad);
        if(sinl2>0){
            l = Math.toDegrees(l2)+l0;
        }else{
            l = Math.toDegrees(2*Math.PI-l2)+l0;
        }
        
        if(l>360){
            l = l-360;
        }
        
        double b = Math.toDegrees(b_rad);
        
        double[] output = {l,b};
        return output;
        
    }    

   /* Converts from SuperGalactic lattitude and longitude (in degrees) to 
   J2000 equatorial coordinates in decimal degrees.*/
    public double[] sgal2Equ(double sgl, double sgb) 
    {
        double sgl0  = -63.55;
        double ra0   = 13.754;
        double obliq = 74.29110005;
        double cosobliq = Math.cos(Math.toRadians(obliq));
        double sinobliq = Math.sin(Math.toRadians(obliq));
        double sgl1     = Math.toRadians(sgl - sgl0);
        
        double sgb_rad = Math.toRadians(sgb);
        
        double dec_rad = Math.asin(Math.sin(sgb_rad)*cosobliq + Math.cos(sgb_rad)*Math.sin(sgl1)*sinobliq);
        double sgl2    = Math.acos(Math.cos(sgb_rad)*Math.cos(sgl1)/Math.cos(dec_rad));
        
        double ra = 0;
        double sinl2 = (Math.cos(sgb_rad)*Math.sin(sgl1)*cosobliq - Math.sin(sgb_rad)*sinobliq)/Math.cos(dec_rad);
        if(sinl2>0){
            ra = Math.toDegrees(sgl2)+ra0;
        }else{
            ra = Math.toDegrees(2*Math.PI-sgl2)+ra0;
        }
        
        if(ra>360){
            ra = ra-360;
        }
        
        double dec = Math.toDegrees(dec_rad);
        
        double[] output = {ra,dec};
        return output;
        
    }    

    /* Converts from J2000 equatorial coordinates in decimal degrees
    to SuperGalactic lattitude and longitude (in degrees).*/
    public double[] equ2Sgal(double ra, double dec) 
    {
        double sgl0  = -63.55;
        double ra0   = 13.754;
        double obliq = 74.29110005;
        double cosobliq = Math.cos(Math.toRadians(obliq));
        double sinobliq = Math.sin(Math.toRadians(obliq));
        double l1    = Math.toRadians(ra - ra0);
        
        double dec_rad = Math.toRadians(dec);
        
        double sgb_rad = Math.asin(Math.sin(dec_rad)*cosobliq - Math.cos(dec_rad)*Math.sin(l1)*sinobliq);
        double l2    = Math.acos(Math.cos(dec_rad)*Math.cos(l1)/Math.cos(sgb_rad));
        
        double sgl = 0;
        double sinl2 = (Math.cos(dec_rad)*Math.sin(l1)*cosobliq + Math.sin(dec_rad)*sinobliq)/Math.cos(sgb_rad);
        if(sinl2>0){
            sgl = Math.toDegrees(l2)+sgl0;
        }else{
            sgl = Math.toDegrees(2*Math.PI-l2)+sgl0;
        }
        
        if(sgl>360){
            sgl = sgl-360;
        }
        
        double sgb = Math.toDegrees(sgb_rad);
        
        double[] output = {sgl,sgb};
        return output;
        
    }    
}