// Horizontal Dial ~ ~ ~ For a given latitude and longitude

// Get latitude, longitude and reference longitude 
lat = input ("Latitude [eg 32.75]: ");
latR = lat * 2 * 3.1416 / 360;
tlatR = tan(latR);
slatR = sin(latR);

lng = input ("Longitude [eg 108.2]: ");
ref = input ("Ref.Longitude [eg 105]: ");
clng = lng - ref;

xset ("font", 1, 8);

// Set the parameters for the graphics area                                 )
xc = 0.65;                                      // x of dial center         ) xc=0.65
yc = 0.1 ;                                      // y of dial center         ) yc = 0.1
ysize = .8 ;                                    // size of the y scale      ) ysize=0.8
xsize = ysize  ;                                // size of the x scale      ) xsize = ysize
                                                //                          )
// This corrects for aspect ratio - must go hand in hand with about values  )
square ( -.2, -.2 , 1.5, 1.5 ) ;                // xmin ymin xmax ymax      ) -.2, -.2, 1.5, 1.5
   
for hr = 7:17  ;                                // 6:19 is normal hour range
   //    However extra code is needed for angles such as for 0600 when dial is 
   //    west of reference longitude, and for angles such as 1800 when the dial
   //    is east of longitude. The purpose of this program was to show the use 
   //    of Scilab's language and techniques.
   //    Simon Wheaton-Smith FRI, MBCS, CITP 
   //    www.illustratingshadows.com             Jan 11, 2008 0815

   // calculate the local hour angle in radians  L H A = local hour angle of the sun
   lha = (15 * hr) - clng    ;                  // local hour angle = hour * 15 degrees
   lhaR = lha * 2 * 3.1416 / 360 ;              // in radians
   tlhaR =  tan(lhaR)         ;                 // tlhaR is tan of hour angle

   // calculate the hour line angle per se       H L A = hour line angle on the plate
   // hrLnAng = atan (  sin(lat) * tan (lha)  )	
   hlaR = atan(slatR * tlhaR)  ;                // hour line angle in radians
   hlaD = ( hlaR * 360) / (2 * 3.1416) ;        // hour line angle in degrees
   
   // determine the end point for early, midday, or late hours
   if ( hlaD < -45 ) then;
      // before 0900
      // assuming x as a fixed end point, what is the y value for the hour
           xe = xc - xsize  ;                   // x end = 1 x size plus dial center x
           ye = yc - (xsize / tan(hlaR));       // y end = x size * formula 
   elseif ( hlaD > 45 ) then;
      // after 1500
      // assuming x as a fixed end point, what is the y value for the hour
      xe = xc + xsize   ;                       // x end = 1 x size plus dial center x
      ye = yc + (xsize / tan(hlaR))  ;          // y end = x size * formula 
    else // 0900 to 1500
      // assuming y as a fixed end point, what is the x value for the hour
      ye = yc + ysize        ;                  // y end = 1 y size plus dial center y
      xe = xc + (ysize * tan(hlaR))  ;          // x end = y size * formula  
   end;
 
   // draw an hour line
   xa = [ xc xe  ];
   ya = [ yc ye  ];
   xarrows(xa,ya);
   
   // name the hour line
   // xstring( 1,  sin(lat)*tan(lha) ,"hr"  )
   xset ("font", 1, 8);
   xstring( xe , ye  , string(hr), 0) ;         // show hour at end of hour line
   xset ("font", 1, 2);
   xstring( xe , ye+0.1  , string(int(hlaD)), 0);
   
end;

// draw a pretty box
xrect  ( xc-xsize, yc+ysize, 2*xsize, ysize) ;  // upper left point, then sizes

// draw a calibration line they must scale to 45 because Scilab aspect ratios 
// leave a lot to be desired.
xa = [ xc       xc+0.2  ] ;                     // center a 45 degree line and use
ya = [ yc-0.2   yc  ] ;                         // x for both x and y to force a 45
xset ("font", 1, 4);
xstring( xc+0.2 , yc-0.1  , string(45), 0) ;    // used for aspect ratio
xarrows(xa,ya);

// state the hour line angles now
printf ("Graphical depiction may or may need aspect ratio adjustment\n");
printf ("Lat: %3.1f  Long.corr was: %3.1f \n",lat, clng);
for hr = 6:18     ;                             // print the hour and hour line angles
    printf ("Hour: %i   Hour Line Angle: %3.1f \n", hr, ( (atan(slatR * (tan((((15 * hr) - clng) * 2 * 3.1416 / 360)))) * 360) / (2 * 3.1416)));
end;
printf ("*** ENTER *** to end program\n");
// END


