'********************************************************
'* NOT FINISHED YET - VARIANT FOR TESTING               *
'* Az-Ht.BAS is a DeltaCad macro for producing a        *
'* Graph with the Azimuth vs. the Height of the Sun     *
'* Created by Valentin Hristov (valhrist@bas.bg) in     *
'* June 2008.                                           *
'* (www.math.bas.bg/complan/valhrist/mystuff.htm)       *
'* I was inspired by Mac Oglesby to use the             *
'* North American Sundial Society DeltaCad programs     *
'* as tutorials (http://sundials.org)                   *
'* and made with DeltaCad (www.deltacad.com)            *
'* different types of sundials.                         *
'********************************************************

'* You can see the picture of my sundial generated by one of my
'* DeltaCad macros at
'* www.flickr.com/photos/Valentin_Hristov/261303801/
'* Click on the button "All sizes" to see a bigger photo with details.
'* I am very grateful to my friends Daniela (www.danyo.net) and
'* Todor (www.todor.org) who converted my drawing into a real art piece!!!

'* There is a Demo version of DeltaCad at www.deltacad.com.
'* Use the menu "Options - Macro - Run..." or the separate "Macro"
'* button - "Edit Macro List", add the file, and "Run Macro".
'* In the dialog box use DECIMAL DEGREES. Negative values indicate
'* South for latitude and West for longitude and central meridian.

'* E N J O Y !!!

Option Explicit ' Force all variables to be declared before they are used. No adhoc variables

dcSetLineParms dcBlack, dcSolid, dcThin

Dim l,p,lon,cm,rot,i,pi,d2r,r2d,splin,dm,dmr,w,ll,ud,lol,cs,j As Double
Dim ha,xe,ye,xb,yb,lc,ham,ham1,hac,hacr As Double
Dim d,h,lcr,slc,dc,clc,lr,sl,cl,dr,sd,cd,ir,si,ci,xx,yy,x,y,hh As Double
Dim rr,sr,cr,dsr,sds,cds,surf,high,feot,flon,bar,ear As Double
Dim x0,y0,z0,x1,y1,z1,x2,y2,z2,x3,y3,z3,x4,y4,z4,shac,chac As Double
Dim r,at,ab,al,ar,per,Time_int,eot_c,lon_c,rrrr,rrr,spring,bb,ee As Double
Dim xnod,ynod,rnod,an,anr,aaa,bbb,ccc,phi,lambda,z,nod,xxx,yyy As Double
Dim coord_c,ovrhd_c,xx1,yy1,xx2,yy2,xlabel,ylabel,sur,abb,att,jj,ang As Double
Dim a,b,lb,la,ba,ea,ovrhd as Double

Dim decl(366),eot(366),spl(732) As Double

Dim button,count,nl,bm(13),ind(37),nom,ind1 As Integer

Dim action,outtext As String

Dim datetext(13),dirs(4) As String

Dim hhv,hhb,hhe,nlcheck As Boolean


dcSetLineParms dcBLACK, dcSOLID, dcTHIN
dcSetCircleParms dcBLACK, dcSOLID, dcTHIN
dcSetTextParms dcPurple, "Tahoma","Bold",0,12,21,0,0

'Establish the 5 standard line thicknesses in thousands of an inch.
dcSetDrawingData dcLineThin,   .003
dcSetDrawingData dcLineNormal, .008
dcSetDrawingData dcLineThick,  .012
dcSetDrawingData dcLineHeavy,  .024
dcSetDrawingData dcLineWide,   .048

'Maximize the window, close any existing drawing without saving, and start a new drawing.
dcSetDrawingWindowMode dcMaximizeWin
dcCloseWithoutSaving
dcNew ""

'dcSetDrawingScale 25.4 'mm

'**************************************

'Start of program
Call Sundial
'**********
Sub Sundial
init_constants
Input_constants_of_sundial
If button=0 Then GoTo cncl
Angles
Main
Height
dcViewAll
cncl:
End Sub
'End of program


'**************************************
'Start of subroutines
'''''''''''''''''''''''''''''''''''''''
Sub Input_constants_of_sundial

Begin Dialog CONSTANTS_INPUT 10,25,205,155, "Input data for the sundial"
 Text 15,0,220,10, "SUN POSITION - AZIMUTH vs HEIGHT (elevation)"
 Text 15,10,250,10, "with longitude and EOT corrections"
 Text 15,20,150,10, "Place"
 TextBox 88,20,99,10, .p
 Text 15,30,180,10, "Latitude(N>0)  Longitude(E>0)  Central meridian(E>0)"
 TextBox 15,40,30,10, .l
 TextBox 65,40,30,10, .lon
 TextBox 115,40,30,10, .cm
 Text 15,50,150,10, "  Use decimal degrees !!!"

 Text 15,70,50,10, "Period"
 OptionGroup .per
 OptionButton 38,60,85,10, "21 Dec - 21 Jun"
 OptionButton 38,70,85,10, "21 Jun - 21 Dec"
 OptionButton 38,80,85,10, "Whole year"

 Text 105,70,80,10, "Time interval"
 OptionGroup .time_int
 OptionButton 150,60,35,10, "15 min"
 OptionButton 150,70,35,10, "30 min"
 OptionButton 150,80,35,10, "60 min"

 Text 15,95,80,10, "EOT correction"
 OptionGroup .eot_c
 OptionButton 65,90,30,10, "Yes"
 OptionButton 65,100,30,10, "No"

 Text 95,95,80,10, "Longitude correction"
 OptionGroup .lon_c
 OptionButton 165,90,30,10, "Yes"
 OptionButton 165,100,30,10, "No"

 Text 15,115,80,10, "Coordinates"
 OptionGroup .coord_c
 OptionButton 65,110,50,10, "Cartesian"
 OptionButton 65,120,75,10, "Polar  ->  Overhead"

 OptionGroup .ovrhd_c
 OptionButton 140,120,25,10, "No"
 OptionButton 165,120,25,10, "Yes"

 OKButton 45,135,37,12
 CANCELButton 115,135,37,12
End Dialog

'Initialize 
Dim prompt As constants_input

'prompt.p = "Sydney - Australia"
'prompt.l = -33.853
'prompt.lon = 151.21
'prompt.cm = 150
prompt.p = "Lozen - Sofia - Bulgaria"
prompt.l = 42.6021
prompt.lon = 23.5011
prompt.cm = 30
prompt.per = 0
prompt.time_int = 0
prompt.eot_c = 0
prompt.lon_c = 0
prompt.coord_c = 1
prompt.ovrhd_c = 0

repeat_until_inputcorrect: 'label to return if input is not correct

button = Dialog(prompt)    'get the input

If (button=0) Then GoTo End_input

If test("Latitude",prompt.l,-90,90) = false Then
  GoTo repeat_until_inputcorrect
End If
If test("Longitude",prompt.lon,-180,180) = false Then
  GoTo repeat_until_inputcorrect
End If
If test("CentralMeridian",prompt.cm,-180,180) = false Then
  GoTo repeat_until_inputcorrect
End If

'Set program variables with input variables, angles in degrees
p = prompt.p
l = CDbl(prompt.l)
lon = CDbl(prompt.lon)
cm = CDbl(prompt.cm)
d = 0
i = 0
rot = 0
per = prompt.per
Time_int=prompt.time_int
eot_c=prompt.eot_c
lon_c=prompt.lon_c
coord_c=prompt.coord_c
ovrhd_c=prompt.ovrhd_c
r=1
al=-r*pi
ar=r*pi
lb=10
la=10
rrrr=1
high=1
If high=0 Then rrr=r-rrrr Else rrr=r+rrrr
at=la
ab=-lb-.3
surf=0

nod=0
rnod=0
an=0
if coord_c=1 then splin=1 else splin=0

If Time_int=0 Then nl=4
If Time_int=1 Then nl=2
If Time_int=2 Then nl=1
If eot_c=0 Then feot=1 Else feot=0
If lon_c=0 Then flon=1 Else flon=0
If at<=ab Or at<0 Or ab>0 Then GoTo repeat_until_inputcorrect
If ar<=al Or al>0 Or ar<0 Then GoTo repeat_until_inputcorrect

ba=-180
ea=180

End_input:

End Sub

'''''''''''''''''''''''''''''''''''''''
Sub init_constants
pi = 4 * Atn(1)
d2r = pi/180
r2d = 180/pi
dm = 23.43954
dmr = dm*d2r

ind( 1)=355
ind( 2)=  1 '1Jan
ind( 3)= 11
ind( 4)= 21
ind( 5)= 32 '1Feb
ind( 6)= 42
ind( 7)= 52
ind( 8)= 60 '1Mar
ind( 9)= 70
ind(10)= 80
ind(11)= 91 '1Apr
ind(12)=101
ind(13)=111
ind(14)=121 '1May
ind(15)=131
ind(16)=141
ind(17)=152 '1Jun
ind(18)=162
ind(19)=172
ind(20)=182 '1Jul
ind(21)=192
ind(22)=202
ind(23)=213 '1Aug
ind(24)=223
ind(25)=233
ind(26)=244 '1Sep
ind(27)=254
ind(28)=264
ind(29)=274 '1Oct
ind(30)=284
ind(31)=294
ind(32)=305 '1Nov
ind(33)=315
ind(34)=325
ind(35)=335 '1Dec
ind(36)=345
ind(37)=355

bm( 1)=  1 '1jan
bm( 2)= 32 '1feb
bm( 3)= 60
bm( 4)= 91
bm( 5)=121
bm( 6)=152
bm( 7)=182
bm( 8)=213
bm( 9)=244
bm(10)=274
bm(11)=305
bm(12)=335 '1dec
bm(13)=366 '1jan

For count=1 To 365
w=.017202792*(count-(cm-15)/360)
decl(count)=.367402-23.275*Cos(w+.178044)-.38506*Cos(2*w+.0687076)-.16046*Cos(3*w+.451301)+.00315469*Cos(4*w+.876643)
Next count
decl(366)=decl(1)

For count=1 To 365
w=.017202792*(count-(cm-15)/360)
eot(count)=.00884207-7.36034*Cos(w+1.49487)+9.91068*Cos(2*w-1.21808)-.306344*Cos(3*w+1.82055)+.204774*Cos(4*w-.875115)
Next count
eot(366)=eot(1)

datetext( 1) = "JAN"
datetext( 2) = "FEB"
datetext( 3) = "MAR"
datetext( 4) = "APR"
datetext( 5) = "MAY"
datetext( 6) = "JUN"
datetext( 7) = "JUL"
datetext( 8) = "AUG"
datetext( 9) = "SEP"
datetext(10) = "OCT"
datetext(11) = "NOV"
datetext(12) = "DEC"
datetext(13) = "JAN"

End Sub

Sub Angles

lc=lon-cm    'longitude correction
lcr=lc*d2r   'longitude correction in radians
slc=Sin(lcr) 'sin(longitude correction)
clc=Cos(lcr) 'cos(longitude correction)

lr=l*d2r     'longitude in radians
sl=Sin(lr)   'sin(latitude)
cl=Cos(lr)   'cos(latitude)

dr=d*d2r     'declination in radians
sd=Sin(dr)   'sin(declination)
cd=Cos(dr)   'cos(declination)

ir=i*d2r     'inclination in radians
si=Sin(ir)   'sin(inclination)
ci=Cos(ir)   'cos(inclination)

rr=rot*d2r   'rotation in radians
sr=Sin(rr)   'sin(rotation)
cr=Cos(rr)   'cos(rotation)

End Sub

Sub Convert

'rotation to horizontal plane
x1=x0
y1=sl*y0+cl*z0
z1=-cl*y0+sl*z0
'rotation for declination of the plane
x2=cd*x1-sd*y1
y2=sd*x1+cd*y1
z2=z1
'rotation for the inclination of the plane
x3=x2
y3=ci*y2+si*z2
z3=-si*y2+ci*z2
'rotation of the box in its plane
x4=cr*x3+sr*y3
y4=-sr*x3+cr*y3
z4=z3

End Sub

'''''''''''''''''''''''''''''''''''
Sub Main

hh=2.5 'Scale for the vertical size

anr=an*d2r  'in radians
xnod=rnod*Sin(anr) 'nodus point
ynod=rnod*Cos(anr)

''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
dcAddLayer "DateLines"
dcSetLineParms dcBlack, dcSolid, dcThin

abb=0
att=0
If per=0 Then
 spring=1
 Drawing
End If
If per=1 Then
 spring=0
 Drawing
End If
If per=2 Then
 spring=1
 Drawing
 spring=0
 Drawing
End If


'Cartesian Coordinate System
if coord_c=0 then

ab=0
at=r*pi/2*hh

dcSetLineParms dcBlack, dcSolid, dcThin
dcSetTextParms dcDarkPurple,"Tahoma","Bold",0,6,21,0,0

for j=0 to 360 step 15
 if l<0 then
  if j>=180 then jj=j-180 else jj=j+180
 else
  jj=j
 end if
 dcCreateLine al+(ar-al)*j/360,ab,al+(ar-al)*j/360,at
 dcCreateText al+(ar-al)*j/360,ab-r*.1,0,CStr(jj)
next j
for j=0 to 360 step 5
 dcCreateLine al+(ar-al)*j/360,ab,al+(ar-al)*j/360,ab-.05
next j
for j=0 to 360
 dcCreateLine al+(ar-al)*j/360,ab,al+(ar-al)*j/360,ab-.025
next j
for j=0 to 90 step 5
 dcCreateLine al,ab+(at-ab)*j/90,ar,ab+(at-ab)*j/90
 dcCreateText al-r*.1,ab+(at-ab)*j/90,0,CStr(j)
 dcCreateText ar+r*.1,ab+(at-ab)*j/90,0,CStr(j)
next j
for j=0 to 90
 dcCreateLine al,ab+(at-ab)*j/90,al-.03,ab+(at-ab)*j/90
 dcCreateLine ar,ab+(at-ab)*j/90,ar+.03,ab+(at-ab)*j/90
next j

dcCreateLine al,0,ar,0
dcCreateLine al,ab,al,at
dcCreateLine al,at,ar,at
dcCreateLine ar,at,ar,ab
dcCreateLine ar,ab,al,ab

bar=ba*d2r
ear=ea*d2r
While bar<-pi
 bar=bar+2*pi
Wend
While bar>pi
 bar=bar-2*pi
Wend
If surf=1 Then bar=-bar
While ear<-pi
 ear=ear+2*pi
Wend
While ear>pi
 ear=ear-2*pi
Wend
If surf=1 Then ear=-ear

dcSetLineParms dcPurple, dcSolid, dcThick
dcCreateLine r*bar,ab,r*bar,at
dcCreateLine r*ear,ab,r*ear,at
dcSetLineParms dcBlack, dcSolid, dcThin

If l>=0 Then
 ll=CStr(l)+" N"
 ud=0
Else
 ll=CStr(-Val(l))+" S"
 ud=180
End If
If lon>=0 Then lol=CStr(lon)+" E" Else lol=CStr(-Val(lon))+" W"
dcSetTextParms dcDarkPurple, "Tahoma","Bold",0,8,21,0,0
dcCreateText (al+ar)/2,ab-0.5,0,p
dcCreateText (al+ar)/2-2.5*r,ab-0.5,0,"Latitude "+ll
dcCreateText (al+ar)/2+2.5*r,ab-0.5,0,"Longitude "+lol
if l<0 then
 dcCreateText al,ab-0.25,0,"S"
 dcCreateText al/2,ab-0.25,0,"W"
 dcCreateText 0,ab-0.25,0,"N"
 dcCreateText ar/2,ab-0.25,0,"E"
 dcCreateText ar,ab-0.25,0,"S"
else
 dcCreateText al,ab-0.25,0,"N"
 dcCreateText al/2,ab-0.25,0,"E"
 dcCreateText 0,ab-0.25,0,"S"
 dcCreateText ar/2,ab-0.25,0,"W"
 dcCreateText ar,ab-0.25,0,"N"
end if
dcCreateText 0,ab-0.35,0,"AZIMUTH"
dcSetTextParms dcDarkPurple,"Tahoma","Bold",90,8,21,0,0
dcCreateText al-.25,(at-ab)/2,0,"HEIGHT"
dcSetTextParms dcDarkPurple,"Tahoma","Bold",0-90,8,21,0,0
dcCreateText ar+.25,(at-ab)/2,0,"HEIGHT"
dcSetTextParms dcDarkPurple,"Tahoma","Bold",0,8,21,0,0
dcCreateText 0,at+.1,0,"SUN POSITION"


If flon=1 Then
 dcCreateText al/3,ab-.65,0,"With Longitude correction"
Else
 dcCreateText al/3,ab-.65,0,"Without Longitude correction"
End If
If feot=1 Then
 dcCreateText ar/3,ab-.65,0,"With EOT correction"
Else
 dcCreateText ar/3,ab-.65,0,"Without EOT correction"
End If


dcSetTextParms dcBlack,"Tahoma","Bold",0,8,21,0,0
dcCreateText (al+ar)/2,ab-.85,0,"Author: Valentin Hristov, Sofia, Bulgaria;       E-mail: valhrist@math.bas.bg"
dcCreateText (al+ar)/2,ab-1,0,"Web page: www.math.bas.bg/complan/valhrist/mystuff.htm"

end if ' Cartesian

dcSetLineParms dcBlack, dcSolid, dcThin

'Polar Coordinate system
if coord_c=1 then

dcCreateLine 1.25*hh*r*pi/2,1.25*hh*r*pi/2,-1.25*hh*r*pi/2,1.25*hh*r*pi/2
dcCreateLine -1.25*hh*r*pi/2,1.25*hh*r*pi/2,-1.25*hh*r*pi/2,-1.25*hh*r*pi/2
dcCreateLine -1.25*hh*r*pi/2,-1.25*hh*r*pi/2,1.25*hh*r*pi/2,-1.25*hh*r*pi/2
dcCreateLine 1.25*hh*r*pi/2,-1.25*hh*r*pi/2,1.25*hh*r*pi/2,1.25*hh*r*pi/2
dcCreateCircle 0,0,1.15*hh*r*pi/2

if ovrhd_c=0 then ovrhd=1 else ovrhd=-1

for i=90 to -269 step -10
 aaa=hh*r*pi/2*cos((i)*d2r)
 bbb=hh*r*pi/2*sin((i)*d2r)
 dcCreateline 0,0,ovrhd*hh*r*pi/2*cos(i*d2r),hh*r*pi/2*sin(i*d2r)
 ''dcSetTextParms dcBlack,"Tahoma","Bold",180+ovrhd*(90-i),10,21,0,0
 dcSetTextParms dcBlack,"Tahoma","Bold",ovrhd*(i+90),10,21,0,0
 dcCreateText 1.05*ovrhd*hh*r*pi/2*cos(i*d2r),1.05*hh*r*pi/2*sin(i*d2r),0,CStr(90-i)
 if i=90 then dcCreateText 1.1*ovrhd*hh*r*pi/2*cos(i*d2r),1.1*hh*r*pi/2*sin(i*d2r),0,"NORTH"
 if i=0 then dcCreateText 1.1*ovrhd*hh*r*pi/2*cos(i*d2r),1.1*hh*r*pi/2*sin(i*d2r),0,"EAST"
 if i=-90 then dcCreateText 1.1*ovrhd*hh*r*pi/2*cos(i*d2r),1.1*hh*r*pi/2*sin(i*d2r),0,"SOUTH"
 if i=-180 then dcCreateText 1.1*ovrhd*hh*r*pi/2*cos(i*d2r),1.1*hh*r*pi/2*sin(i*d2r),0,"WEST"
next i

for i=10 to 90 step 10
 dcCreateCircle 0,0,i*hh*r*d2r
 dcSetTextParms dcBlack,"Tahoma","Bold",45,10,21,0,0
 dcCreateText i*hh*r*d2r*cos(-pi/4),i*hh*r*d2r*sin(-pi/4),0,CStr(90-i)
 dcSetTextParms dcBlack,"Tahoma","Bold",135,10,21,0,0
 dcCreateText i*hh*r*d2r*cos(pi/4),i*hh*r*d2r*sin(pi/4),0,CStr(90-i)
 dcSetTextParms dcBlack,"Tahoma","Bold",215,10,21,0,0
 dcCreateText i*hh*r*d2r*cos(3*pi/4),i*hh*r*d2r*sin(3*pi/4),0,CStr(90-i)
 dcSetTextParms dcBlack,"Tahoma","Bold",315,10,21,0,0
 dcCreateText i*hh*r*d2r*cos(5*pi/4),i*hh*r*d2r*sin(5*pi/4),0,CStr(90-i)
next i

dcSetTextParms dcBlack,"Tahoma","Bold",0,12,21,0,0
dcCreateCText 0,-1.2*hh*r*pi/2,-1.2*hh*r*pi/2,"SUN POSITION"
dcSetTextParms dcBlack,"Tahoma","Bold",180,12,21,0,0
dcCreateCText 0,1.2*hh*r*pi/2,-1.2*hh*r*pi/2,p
If l>=0 Then
 ll=CStr(l)+" N"
 ud=0
Else
 ll=CStr(-Val(l))+" S"
 ud=180
End If
dcSetTextParms dcBlack,"Tahoma","Bold",135,12,21,0,0
dcCreateCText 1.2*hh*r*pi/2/sqr(2),1.2*hh*r*pi/2/sqr(2),-1.2*hh*r*pi/2,"Latitude   "+ll
If lon>=0 Then lol=CStr(lon)+" E" Else lol=CStr(-Val(lon))+" W"
dcSetTextParms dcBlack,"Tahoma","Bold",-135,12,21,0,0
dcCreateCText -1.2*hh*r*pi/2/sqr(2),1.2*hh*r*pi/2/sqr(2),-1.2*hh*r*pi/2,"Longitude   "+lol
dcSetTextParms dcBlack,"Tahoma","Bold",-90,12,21,0,0
dcCreateText -1.2*hh*r*pi/2,0,0,"Valentin Hristov   E-mail: valhrist@math.bas.bg"
dcSetTextParms dcBlack,"Tahoma","Bold",90,12,21,0,0
dcCreateText 1.2*hh*r*pi/2,0,0,"Web page: www.math.bas.bg/complan/valhrist/mystuff.htm"

dcSetTextParms dcBlack,"Tahoma","Bold",-45,12,21,0,0
If flon=1 Then
 dcCreateCText -1.2*hh*r*pi/2/sqr(2),-1.2*hh*r*pi/2/sqr(2),-1.2*hh*r*pi/2,"With Longitude correction"
Else
 dcCreateCText -1.2*hh*r*pi/2/sqr(2),-1.2*hh*r*pi/2/sqr(2),-1.2*hh*r*pi/2,"Without Longitude correction"
End If
dcSetTextParms dcBlack,"Tahoma","Bold",45,12,21,0,0
If feot=1 Then
 dcCreateCText 1.2*hh*r*pi/2/sqr(2),-1.2*hh*r*pi/2/sqr(2),-1.2*hh*r*pi/2,"With EOT correction"
Else
 dcCreateCText 1.2*hh*r*pi/2/sqr(2),-1.2*hh*r*pi/2/sqr(2),-1.2*hh*r*pi/2,"Without EOT correction"
End If

end if ' Polar

End Sub 'Main

Sub Drawing
' hour lines
cs=0                         'counter for spline
For ha=0 To 23.999 Step 1/nl
 ' civil and daylight savings labels
 ham=ha-1
 ham=ham-Int(ham/12)*12+1 'ha (12)
 ham1=ham
 ham1=ham1-Int(ham1/12)*12+1 'ha+1 (12)

 hhv=False
 hhb=True
 hhe=False
 If spring=1 Then nom=183 Else nom=184
 For ind1=1 To nom
  If spring=1 Then
   If ind1<12 Then count=354+ind1 Else count=ind1-11
  Else
   count=171+ind1
  End If
  hac=ha*15+flon*lc-feot*eot(count)/4
  hacr=hac*d2r
  shac=Sin(hacr)
  chac=Cos(hacr)
  dsr=decl(count)*d2r
  sds=Sin(-dsr)
  cds=Cos(-dsr)
  x0=-shac*cds
  y0=-chac*cds
  z0=sds

  convert

If nod=1 Then
 Nodus
Else
If x4>0 Then
 xx=r*(pi/2-arcsin(y4/Sqr(x4*x4+y4*y4)))
Else
 xx=r*(-pi/2-arcsin(-y4/Sqr(x4*x4+y4*y4)))
End If

if l<0 and coord_c=0 then
 if xx<0 then
  xx=xx+r*pi
 else
  xx=xx-r*pi
 end if
end if

 yy=Atn(-z4/Sqr(x4*x4+y4*y4))*r*hh
End If
If surf=1 Then xx=-xx

if coord_c=1 then 'polar
 xxx=(hh*pi/2-yy)*cos(xx+pi/2)
 yyy=(hh*pi/2-yy)*sin(xx+pi/2)
 yy=yyy
 if ovrhd_c=0 then
  xx=0-xxx
 else
  xx=xxx
 end if
end if

'  cylinder

   if z1<0 then
'  If (z1<0 _
'      And xx>al And xx<ar _
'      And yy>ab And yy<at) Then 'illuminated rectangle
   cs=cs+1
   spl(2*cs-1)=xx
   spl(2*cs)=yy
   If yy<abb Then abb=yy
   If yy>att Then att=yy

   xe=xx
   ye=yy
   If hhb=True Then
    xb=xx
    yb=yy
    hhb=False
    hhv=True
   End If
  Else
   If cs>2 Then
    If spring=1 Then
     If ha=Int(ha) Then
      dcSetLineParms dcBlue, dcSOLID, dcTHICK
      dcSetSplineParms dcBlue, dcSOLID, dcTHICK
     Else
      dcSetLineParms dcBlue, dcSOLID, dcTHIN
      dcSetSplineParms dcBlue, dcSOLID, dcTHIN
     End If
    Else
     If ha=Int(ha) Then
      dcSetLineParms dcRed, dcSOLID, dcTHICK
      dcSetSplineParms dcRed, dcSOLID, dcTHICK
     Else
      dcSetLineParms dcRed, dcSOLID, dcTHIN
      dcSetSplineParms dcRed, dcSOLID, dcTHIN
     End If
    End If

    Draw

    cs=0
   Else
    If cs=2 Then
     dcCreateLine spl(1),spl(2),spl(3),spl(4)
    End If
    cs=0
   End If
  End If ' illuminated
 Next ind1

 hhe=False
 If cs>2 Then                       'if the hour line finishes at 365
  If spring=1 Then
   If ha=Int(ha) Then
    dcSetLineParms dcBlue, dcSOLID, dcTHICK
    dcSetSplineParms dcBlue, dcSOLID, dcTHICK
   Else
    dcSetLineParms dcBlue, dcSOLID, dcTHIN
    dcSetSplineParms dcBlue, dcSOLID, dcTHIN
   End If
  Else
   If ha=Int(ha) Then
    dcSetLineParms dcRed, dcSOLID, dcTHICK
    dcSetSplineParms dcRed, dcSOLID, dcTHICK
   Else
    dcSetLineParms dcRed, dcSOLID, dcTHIN
    dcSetSplineParms dcRed, dcSOLID, dcTHIN
   End If
  End If
  Draw
 Else
  If cs=2 Then
   dcCreateLine spl(1),spl(2),spl(3),spl(4)
  End If
 End If

 'Cartesian
 If coord_c=0 then

 If (spring=1 And l>=0) Or (spring=0 And l<0) Then
  If ha=Int(ha) And hhe=False And hhv=True Then
   dcSetTextParms dcBlue, "Tahoma","Bold",0,8,21,0,0
   dcCreateText xe,ye+.1,0,CStr(ham)
   dcSetTextParms dcRed, "Tahoma","Bold",0,8,21,0,0
   dcCreateText xe,ye+.2,0,CStr(ham1)
  End If
 Else
  If ha=Int(ha) And hhe=False And hhv=True Then
   dcSetTextParms dcBlue, "Tahoma","Bold",0,8,21,0,0
   dcCreateText xb,yb+.1,0,CStr(ham)
   dcSetTextParms dcRed, "Tahoma","Bold",0,8,21,0,0
   dcCreateText xb,yb+.2,0,CStr(ham1)
  End If
 End If
 End If

 'Polar
 If coord_c=1 then

 If spring=1 Or spring=2 Then
  If ha=Int(ha) And hhe=False And hhv=True Then
   dcSetTextParms dcBlue, "Tahoma","Bold",0,8,21,0,0
   dcCreateText xb,yb+.1,0,CStr(ham)
   dcSetTextParms dcRed, "Tahoma","Bold",180,8,21,0,0
   dcCreateText xe,ye-.1,0,CStr(ham1)
  End If
 Else
  If ha=Int(ha) And hhe=False And hhv=True Then
   dcSetTextParms dcBlue, "Tahoma","Bold",0,8,21,0,0
   dcCreateText xe,ye+.1,0,CStr(ham)
   dcSetTextParms dcRed, "Tahoma","Bold",180,8,21,0,0
   dcCreateText xb,yb-.1,0,CStr(ham1)
  End If
 End If
 End If

 cs=0

Next ha

''''''''''''''''''''''''''''''''''''''''''''''''
dcSetCurrentLayer "DateLines"
''''''''''''''''''''''''''''''''''''''''''''''''
'Date lines
cs=0                         'counter for spline
If spring=1 Then
 bb=1
 ee=19
Else
 bb=19
 ee=37
End If

For j=bb To ee
 xlabel=hh*r*pi
 ylabel=0
 count=ind(j)
 hhv=False
 hhb=True
 hhe=False
 For ha=0 To 23.999 Step .1
  hac=ha*15+flon*lc-feot*eot(count)/4
  hacr=hac*d2r
  shac=Sin(hacr)
  chac=Cos(hacr)
  dsr=decl(count)*d2r
  sds=Sin(-dsr)
  cds=Cos(-dsr)
  x0=-shac*cds
  y0=-chac*cds
  z0=sds

  convert

If nod=1 Then
 Nodus
Else
 If x4>0 Then
  xx=r*(pi/2-arcsin(y4/Sqr(x4*x4+y4*y4)))
 Else
  xx=r*(-pi/2-arcsin(-y4/Sqr(x4*x4+y4*y4)))
 End If

if l<0 and coord_c=0 then
 if xx<0 then
  xx=xx+r*pi
 else
  xx=xx-r*pi
 end if
end if

 yy=Atn(-z4/Sqr(x4*x4+y4*y4))*r*hh
End If

If surf=1 Then xx=-xx

if coord_c=1 then 'polar
 xxx=(hh*pi/2-yy)*cos(xx+pi/2)
 yyy=(hh*pi/2-yy)*sin(xx+pi/2)
 yy=yyy
 if ovrhd_c=0 then
  xx=0-xxx
 else
  xx=xxx
 end if
end if

   if z1<0 then
'  If (z1<0 _
'      And xx>al And xx<ar _
'      And yy>ab And yy<at) Then 'illuminated rectangle
   cs=cs+1
   If abs(xx)<abs(xlabel) Then
    xlabel=xx
    ylabel=yy
   End If
   spl(2*cs-1)=xx
   spl(2*cs)=yy
   If yy<abb Then abb=yy
   If yy>att Then att=yy

   xe=xx
   ye=yy
   If hhb=True Then
    xb=xx
    yb=yy
    hhb=False
    hhv=True
   End If
  Else
   If cs>2 Then
    If j+1=Int((j+1)/3)*3 Then
     If spring=1 Then
      dcSetLineParms dcBlue, dcSOLID, dcTHICK
      dcSetSplineParms dcBlue, dcSOLID, dcTHICK
     Else
      dcSetLineParms dcRed, dcSOLID, dcTHICK
      dcSetSplineParms dcRed, dcSOLID, dcTHICK
     End If
    Else
     If spring=1 Then
      dcSetLineParms dcBlue, dcSOLID, dcTHIN
      dcSetSplineParms dcBlue, dcSOLID, dcTHIN
     Else
      dcSetLineParms dcRed, dcSOLID, dcTHIN
      dcSetSplineParms dcRed, dcSOLID, dcTHIN
     End If
    End If

   Draw

    cs=0
   Else
    If cs=2 Then
     dcCreateLine spl(1),spl(2),spl(3),spl(4)
    End If
    cs=0
   End If
  End If
 Next ha



 hhe=False
 If cs>2 Then
    If j+1=Int((j+1)/3)*3 Then
     If spring=1 Then
      dcSetLineParms dcBlue, dcSOLID, dcTHICK
      dcSetSplineParms dcBlue, dcSOLID, dcTHICK
     Else
      dcSetLineParms dcRed, dcSOLID, dcTHICK
      dcSetSplineParms dcRed, dcSOLID, dcTHICK
     End If
    Else
     If spring=1 Then
      dcSetLineParms dcBlue, dcSOLID, dcTHIN
      dcSetSplineParms dcBlue, dcSOLID, dcTHIN
     Else
      dcSetLineParms dcRed, dcSOLID, dcTHIN
      dcSetSplineParms dcRed, dcSOLID, dcTHIN
     End If
    End If
  Draw
 Else
  If cs=2 Then
   dcCreateLine spl(1),spl(2),spl(3),spl(4)
  End If
 End If
 cs=0
  If j+1=Int((j+1)/3)*3 Then
   If coord_c=1 then
    If l>=0 then
     If spring=1 then
      dcSetTextParms dcBlue, "Tahoma","Bold",180,10,21,0,0
     Else
      dcSetTextParms dcRed, "Tahoma","Bold",180,10,21,0,0
     End If
    Else
     If spring=1 then
      dcSetTextParms dcRed, "Tahoma","Bold",180,10,21,0,0
     Else
      dcSetTextParms dcBlue, "Tahoma","Bold",180,10,21,0,0
      End If
     End If
   End If
   'if ylabel>0 then
   dcCreateText 0,ylabel,0,"1 "+datetext((j+1)/3)
  End If

Next j

dcSetCurrentLayer "default"

End Sub 'Drawing

'''''''''''''''''''''''''''''''''''''''
Sub Height

dcCreateLine -1.25*hh*r*pi/2,-2.25*hh*r*pi/2,1.25*hh*r*pi/2,-2.25*hh*r*pi/2
dcCreateLine -1.25*hh*r*pi/2,-1.25*hh*r*pi/2,-1.25*hh*r*pi/2,-2.25*hh*r*pi/2
dcCreateLine 1.25*hh*r*pi/2,-1.25*hh*r*pi/2,1.25*hh*r*pi/2,-2.25*hh*r*pi/2
dcCreateLine hh*r*pi/2,-1.25*hh*r*pi/2,hh*r*pi/2,-2.25*hh*r*pi/2
'dcCreateLine -hh*r*pi/2,-1.25*hh*r*pi/2,-hh*r*pi/2,-2.25*hh*r*pi/2
'dcCreateLine 0,-1.25*hh*r*pi/2,0,-2.25*hh*r*pi/2

for i=1 to 89
 dcCreateLine _
 hh*r*pi/2-hh*r*pi/2*cos(i*d2r),-2.25*hh*r*pi/2+hh*r*pi/2*sin(i*d2r), _
 hh*r*pi/2-.98*hh*r*pi/2*cos(i*d2r),-2.25*hh*r*pi/2+.98*hh*r*pi/2*sin(i*d2r)
 dcCreateLine _
 -.25*hh*r*pi/2-i*hh*r*pi/2/90,-1.25*hh*r*pi/2, _
 -.25*hh*r*pi/2-i*hh*r*pi/2/90,-1.02*1.25*hh*r*pi/2
next i
for i=0 to 85 step 5
 dcCreateLine _
 hh*r*pi/2-hh*r*pi/2*cos(i*d2r),-2.25*hh*r*pi/2+hh*r*pi/2*sin(i*d2r), _
 hh*r*pi/2-.96*hh*r*pi/2*cos(i*d2r),-2.25*hh*r*pi/2+.96*hh*r*pi/2*sin(i*d2r)
 if i<>0 then
  dcSetTextParms dcBlack,"Tahoma","Bold",180-i,12,21,0,0
  dcCreateText hh*r*pi/2-.92*hh*r*pi/2*cos(i*d2r), _
  -2.25*hh*r*pi/2+.92*hh*r*pi/2*sin(i*d2r),0,CStr(i)
 end if
 dcCreateLine _
 -.25*hh*r*pi/2-i*hh*r*pi/2/90,-1.25*hh*r*pi/2, _
 -.25*hh*r*pi/2-i*hh*r*pi/2/90,-1.04*1.25*hh*r*pi/2
 dcSetTextParms dcBlack,"Tahoma","Bold",180,8,21,0,0
 dcCreateText -.25*hh*r*pi/2-i*hh*r*pi/2/90, _
 -1.06*1.25*hh*r*pi/2,0,CStr(i)
next i

dcSetTextParms dcBlack,"Tahoma","Bold",180,12,0,0,0
dcCreateText -hh*r*pi/16,-2.22*hh*r*pi/2,3.2,"Use the protractor"+ _
 " to find the HEIGHT OF THE SUN."
dcCreateText -hh*r*pi/16,-2.12*hh*r*pi/2,3.2,"Then the main square can"+ _
 " be used as SUNDIAL and SUN COMPASS:"
dcCreateText -hh*r*pi/16,-2.02*hh*r*pi/2,3.2,"1) The TIME can be found"+ _
 " at the intersection point of the circumference with the"+ _
 " determined HEIGHT and the DATE LINE. Use the height ruler to"+ _
 " measure the radius."
dcCreateText -hh*r*pi/16,-1.8*hh*r*pi/2,3.2,"2) If the vector from the"+ _
 " centre to the determined point is directed towards the sun, then"+ _
 " North-South and East-West are correctly shown if the option"+ _
 " <Overhead> is NOT used."
dcCreateText -hh*r*pi/16,-1.57*hh*r*pi/2,3.2,"The <Overhead> option"+ _
 " gives a realistic (not reversed) view of the sun path if you keep"+ _
 " the polar diagram vertically with the desired direction down."
dcCreateText .95*hh*r*pi/2,-2.05*hh*r*pi/2,0,"TOWARDS THE SUN"
dcCreateLine .95*hh*r*pi/2,-2.15*hh*r*pi/2,.5*hh*r*pi/2,-2.15*hh*r*pi/2
dcCreateLine .95*hh*r*pi/2,-2.15*hh*r*pi/2,.85*hh*r*pi/2,-2.1*hh*r*pi/2
dcCreateLine .95*hh*r*pi/2,-2.15*hh*r*pi/2,.85*hh*r*pi/2,-2.2*hh*r*pi/2

dcSetTextParms dcBlack,"Courier New","Normal",180,10,21,0,0
dcCreateText .4*hh*r*pi/2,-1.28*hh*r*pi/2,0,"cut"
dcSetTextParms dcBlack,"Courier New","Normal",90,10,21,0,0
dcCreateText -1.15*hh*r*pi/2,-1.75*hh*r*pi/2,2.25,"Keep vertivally"+ _
 " with this vertex at the centre of the square and rotate."
dcCreateLine -1.25*hh*r*pi/2,-1.25*hh*r*pi/2,-1.15*hh*r*pi/2,-1.45*hh*r*pi/2
dcCreateText 1.03*hh*r*pi/2,-1.75*hh*r*pi/2,0,"fold perpendicularly"
dcCreateText 1.15*hh*r*pi/2,-1.8*hh*r*pi/2,0,"shadow casting point"
dcCreateLine 1.15*hh*r*pi/2,-2.1*hh*r*pi/2,hh*r*pi/2,-2.25*hh*r*pi/2
dcCreateLine 1.1*hh*r*pi/2,-2.2*hh*r*pi/2,hh*r*pi/2,-2.25*hh*r*pi/2
dcCreateLine 1.05*hh*r*pi/2,-2.15*hh*r*pi/2,hh*r*pi/2,-2.25*hh*r*pi/2

End Sub 'Height


'''''''''''''''''''''''''''''''''''''''
Function arcsin(ByVal x) As Double
If Abs(x) > 0.999999999999 Then x = sgn(x)*0.999999999999
  arcsin = Atn(x/Sqr(1-x*x))
End Function

Function arccos(ByVal x) As Double
  arccos = pi/2-arcsin(x)
End Function

Function test(varname,x,minval,maxval) As boolean
If IsNumeric(x) = false Then
test = false
outtext = varname & " must be numeric"
MsgBox outtext
exit Function
End If
If x < minval Or x > maxval Then
outtext = varname & " must be between " & chr$(13) & minval & "  and  " & maxval
MsgBox outtext
exit Function
End If
test = true
End Function

Sub Draw

'splin=1  'does not work properly if lat > 66.5
'splin=0  'draws only the short segments

If splin=0 Then
 For jj=1 To cs-1
  xx1=spl(2*jj-1)
  yy1=spl(2*jj)
  xx2=spl(2*jj+1)
  yy2=spl(2*jj+2)
  If Abs(xx2-xx1)<.2 Then     'only short segments
   dcCreateLine xx1,yy1,xx2,yy2
  End If
 Next jj
Else
 dcCreateSpline spl(1),cs,False
End If
End Sub 'Draw

Sub Nodus
'Variant "azimuth like" angle

aaa=x4*x4+y4*y4        'coefficients of the quadratic equation
bbb=xnod*x4+ynod*y4
ccc=xnod*xnod+ynod*ynod-r*r

'if z4<0 then
lambda=(-bbb+Sqr(bbb*bbb-aaa*ccc))/aaa 'solution of the quadratic equation
'else
'lambda=(-bbb-sqr(bbb*bbb-aaa*ccc))/aaa 'solution of the quadratic equation
'end if
x=xnod+lambda*x4   'intersection of the beam with the cylinder
y=ynod+lambda*y4
z=lambda*z4

'cases
If y=0 Then
 If x>0 Then xx=r*pi/2
 If x<0 Then xx=-r*pi/2
Else
 phi=Atn(x/y)
 If y<0 Then
  If x=0 Then phi=r*pi
  If x>0 Then phi=phi+pi
  If x<0 Then phi=phi-pi
 End If
 xx=r*phi     'flat surface of the cylinder
End If

yy=z

End Sub 'Nodus