rem Pippes ........... Rev 3.2
      rem A J Tooth // March 2006
      on error if (err=17) then quit
      
      himem=lomem + 30000000
      
      rem Setup
      proc_setup
      
      rem Set up rocks
      proc_pipes(1)
      
      a$=""
      *REFRESH OFF
      repeat
        sys "SetStretchBltMode", @memhdc%, 3
        command$="MDISPLAY "+str$~pic%+str$(-xscreen%)+","+str$(-yscreen%)+","+str$(2*xscreen%)+","+str$(2*yscreen%)
        oscli command$
        
        for a%=1 to num%
          rem Display the pipes
          proc_disp(pipe{(a%)})
        next a%
        
        a$=inkey$(10)
        rem Reset pipes
        proc_pipes(0)
        cls
      until a$<>""
      quit
      end
      rem End of Program ==========================================
      rem =========================================================
      
      rem Setup
      def  proc_setup
      *FLOAT 64
      
      rem Go to fullscreen
      proc_fullscreen
      
      A=-45.0 : Dis=2500.0
      
      dim pic% 2359350
      command$=" LOAD "+chr$(34) + @dir$ + "BackPic.bmp" +chr$(34)+" "+str$~pic%
      oscli command$
      
      dim R(3) : rem Radius vector
      
      dim drw% 500, calc% 1500, ftmp% 7, itmp% 3, f50% 7, rgb% 3
      dim Red% 0, Gre% 0, Blu% 0
      
      |f50%=50.0
      
      rem Dual-pass assembly, in case of labels
      for pass&=0 to 2 step 2
        proc_drw(pass&)
        proc_calc(pass&)
      next pass&
      cls
      
      dim light{(2)Ls(3),cl&(3)}
      
      rem Light source vectors
      light{(1)}.Ls(1)=1.0
      light{(1)}.Ls(2)=-1.0
      light{(1)}.Ls(3)=0.0
      light{(2)}.Ls(1)=-1.0
      light{(2)}.Ls(2)=0.0
      light{(2)}.Ls(3)=1.0
      
      rem Light source hues
      light{(1)}.cl&(1)=50
      light{(1)}.cl&(2)=50
      light{(1)}.cl&(3)=50
      light{(2)}.cl&(1)=100
      light{(2)}.cl&(2)=100
      light{(2)}.cl&(3)=0
      
      rem Normalise
      for j&=1 to 2
        S=0.0
        for i&=1 to 3
          S+=(light{(j&)}.Ls(i&))^2
        next i&
        S=sqr(S)
        
        for i&=1 to 3
          light{(j&)}.Ls(i&)=(light{(j&)}.Ls(i&))/S
        next i&
      next j&
      
      num%=25
      origin xscreen%,yscreen%
      endproc
      rem =========================================================
      
      rem Set up use of Full Screen
      def proc_fullscreen
      
      sys "GetSystemMetrics", 0 to xscreen%
      sys "GetSystemMetrics", 1 to yscreen%
      sys "SetWindowLong",@hwnd%,-16,&16000000
      sys "SetWindowPos",@hwnd%,-1,0,0,xscreen%,yscreen%,0
      vdu 23,22,xscreen%;yscreen%;8,16,16,1 : rem Set fullscreen mode
      mouse off : off : rem Turns off the Mouse Pointer and the Cursor
      endproc
      rem -------------------------------------------------------------
      
      rem Norm a vector
      def proc_norm(return Vec())
      local M
      
      M=sqr(Vec(1)*Vec(1)+Vec(2)*Vec(2)+Vec(3)*Vec(3))
      Vec()/=M
      endproc
      rem -------------------------------------------------------------
      
      rem Display the point
      def proc_disp(rr{})
      
      Xs=rr.pst(1) : Zs=rr.pst(2) : Ys=rr.pst(3)
      Xe=rr.pen(1) : Ze=rr.pen(2) : Ye=rr.pen(3)
      
      rem Direction vector for this pipe
      Xd=Xe-Xs : Zd=Ze-Zs : Yd=Ye-Ys
      if Zd>0.0 then Zd=-Zd : Xd=-Xd : Yd=-Yd
      Md=sqr(Xd*Xd + Zd*Zd + Yd*Yd)
      
      rem 'x' and 'y' unit Normals to pipe direction vector
      if abs(Xd)>0.0 then Nx=-1/Xd else Nx=1000.0
      if abs(Zd)>0.0 then Nz=1/Zd else Nz=1000.0
      if abs(Yd)>0.0 then Ny=-1/Yd else Ny=1000.0
      Mx=sqr(Nx*Nx + Nz*Nz) : My=sqr(Ny*Ny + Nz*Nz)
      if Mx=0.0 then Mx=0.001
      if My=0.0 then My=0.001
      Nxx=Nx/Mx : Nxz=Nz/Mx : Nxy=0.0   : rem 'x' normal
      Nyx=0.0   : Nyz=Nz/My : Nyy=Ny/My : rem 'y' normal
      
      rem Draw for each angle
      rem Double-pass
      proc_angle(1)
      proc_angle(0)
      endproc
      rem ===================================================================
      
      rem Set up pipes
      def proc_pipes(Ctr&)
      local a%,j&,i&,Mx%
      
      if Ctr&=1 then
        dim pipe{(num%)pst(3),pen(3),col&(3)}
      endif
      
      for a%=1 to num%
        for j&=1 to 3
          Mx%=50
          rem Position vector - start
          pipe{(a%)}.pst(j&)=1.0*(rnd(2*Mx%)-Mx%)
          
          rem Position vector - end
          pipe{(a%)}.pen(j&)=1.0*(rnd(2*Mx%)-Mx%)
        next j&
        
        repeat
          Zdif=abs(pipe{(a%)}.pen(2) - pipe{(a%)}.pst(2))
          if Zdif<10.0 then pipe{(a%)}.pen(2)=1.0*(rnd(2*Mx%)-Mx%)
        until Zdif>=10.0
        
        rem Colour hue
        for i&=1 to 3
          pipe{(a%)}.col&(i&)=rnd(100)
        next i&
        
      next a%
      endproc
      rem ===================================================================
      
      rem Draw for each angle
      def proc_angle(Ct&)
      local th%,ang,cs,sn,csA,csB
      
      for th%=0 to 360
        ang=th%*2*pi/360
        cs=cos(ang) : sn=sin(ang)
        
        rem Vector to the pipe surface
        R(1)=1.0*cs*Nxx
        R(2)=1.0*cs*Nxz + sn*Nyz
        R(3)=1.0*sn*Nyy
        
        proc_norm(R())
        
        rem Shading
        csA=R(1)*light{(1)}.Ls(1) + R(2)*light{(1)}.Ls(2) + R(3)*light{(1)}.Ls(3)
        csB=R(1)*light{(2)}.Ls(1) + R(2)*light{(2)}.Ls(2) + R(3)*light{(2)}.Ls(3)
        
        ddr&=-int(light{(1)}.cl&(1)*csA*(csA>0.0) + light{(2)}.cl&(1)*csB*(csB>0.0))
        ddg&=-int(light{(1)}.cl&(2)*csA*(csA>0.0) + light{(2)}.cl&(2)*csB*(csB>0.0))
        ddb&=-int(light{(1)}.cl&(3)*csA*(csA>0.0) + light{(2)}.cl&(3)*csB*(csB>0.0))
        
        rs%=rr.col&(1) + ddr& : if rs%>255 then rs%=255
        rf&=rs%
        gs%=rr.col&(2) + ddg& : if gs%>255 then gs%=255
        gf&=gs%
        bs%=rr.col&(3) + ddb& : if bs%>255 then bs%=255
        bf&=bs%
        
        colour 10,rf&,gf&,bf& : gcol 10
        
        rem Set plot points depending on pipe-length
        Lm%=20*int(Md)
        
        sys "GetCurrentProcess" to hprocess%
        sys "SetPriorityClass", hprocess%, &100
        rem Draw the pipe
        if (Ct&=1 and R(2)<0.0) then
          call calc%
        endif
        if Ct&=0 then
          call calc%
        endif
        if th%mod3=0 then *REFRESH
        sys "GetCurrentProcess" to hprocess%
        sys "SetPriorityClass", hprocess%, &20
      next th%
      endproc
      rem ===================================================================
      
      rem Assembly Routine for Plotting
      def proc_drw(opt&)
      P%=drw%
      
      [opt opt&
      
      mov eax,[^yscreen%]
      sub eax,[^ym%]
      shr eax,1
      push eax
      mov eax,[^xscreen%]
      add eax,[^xm%]
      shr eax,1
      push eax
      push @memhdc%
      call "GetPixel"
      mov [rgb%],eax
      mov ecx,0
      mov cl,al
      shr eax,8
      mov edx,0
      mov dl,ah
      sub al,dl
      sub cl,dl
      add al,cl
      cmp al,0
      ja near over
      
      mov bl,[^rf&]
      shr bl,2
      mov dl,bl
      add bl,dl
      add bl,dl
      mov cl,rgb%[2]
      shr cl,2
      add bl,cl
      mov [Red%],bl
      
      mov bl,[^gf&]
      shr bl,2
      mov dl,bl
      add bl,dl
      add bl,dl
      mov cl,rgb%[1]
      shr cl,2
      add bl,cl
      mov [Gre%],bl
      
      mov bl,[^bf&]
      shr bl,2
      mov dl,bl
      add bl,dl
      add bl,dl
      mov cl,[rgb%]
      shr cl,2
      add bl,cl
      mov [Blu%],bl
      
      mov al,19               ;Changes the rgb vlue f0r colr 10
      call "oswrch"
      mov al,10
      call "oswrch"
      mov al,16
      call "oswrch"
      mov al,[Red%]
      call "oswrch"
      mov al,[Gre%]
      call "oswrch"
      mov al,[Blu%]
      call "oswrch"
      
      mov al,18              ;Calls gc0l routine
      call "oswrch"
      mov al,0
      call "oswrch"
      mov al,10
      call "oswrch"
      
      mov al,25
      call "oswrch"
      mov al,69              ;Calls pl0t 69 routine
      call "oswrch"
      mov ebx,[^xm%]
      mov al,bl
      call "oswrch"
      mov al,bh
      call "oswrch"
      mov ebx,[^ym%]
      mov al,bl
      call "oswrch"
      mov al,bh
      call "oswrch"
      
      mov al,25
      call "oswrch"
      mov al,5          ;Calls pl0t 5 routine
      call "oswrch"
      mov ebx,[^xm%]
      mov al,bl
      call "oswrch"
      mov al,bh
      call "oswrch"
      mov ebx,[^ym%]
      mov al,bl
      call "oswrch"
      mov al,bh
      call "oswrch"
      
      .over ret
      ]
      endproc
      
      rem ===================================================================
      
      rem Assembly Routine for Calculation
      def proc_calc(opt&)
      P%=calc%
      
      [opt opt&
      
      mov edx,0
      mov [^c%],edx
      
      .loop
      finit
      fild dword [^c%]
      fild dword [^Lm%]
      fdivp st1,st0
      fstp qword [^ff]
      
      fld qword [^ff]
      fld1
      fsubrp st1,st0
      fld qword [^Zs]
      fmulp st1,st0
      fld qword [^ff]
      fld qword [^Ze]
      fmulp st1,st0
      faddp st1,st0
      fstp qword [^Zmm]
      
      fld qword [^A]
      fld1
      fsubrp st1,st0
      fld qword [^Zmm]
      fld qword [^Dis]
      faddp st1,st0
      fld qword [^A]
      fsubp st1,st0
      fdivp st1,st0
      fstp qword [^Fac]
      
      fld qword [^Fac]
      fld st0
      fld qword [^Xs]
      fmulp st1,st0
      fstp qword [^xs]
      fld st0
      fld qword [^Ys]
      fmulp st1,st0
      fstp qword [^ys]
      fld st0
      fld qword [^Xe]
      fmulp st1,st0
      fstp qword [^xe]
      fld qword [^Ye]
      fmulp st1,st0
      fstp qword [^ye]
      
      fld qword [f50%]
      fld qword [^Zmm]
      fsubp st1,st0
      fistp dword [^rdl%]
      
      fld qword [^Xd]             ;Calculate xm%
      fld st0
      fmulp st1,st0
      fld qword [^Md]
      fld st0
      fmulp st1,st0
      fsubrp st1,st0
      fsqrt
      fld qword [^Md]
      fdivp st1,st0
      
      fld qword [^R(1)]
      fmulp st1,st0
      fild dword [^rdl%]
      fmulp st1,st0
      fstp qword [ftmp%]
      
      fld qword [^ff]
      fld1
      fsubrp st1,st0
      fld qword [^xs]
      fmulp st1,st0
      fld qword [^ff]
      fld qword [^xe]
      fmulp st1,st0
      faddp st1,st0
      fild dword [^xscreen%]
      fmulp st1,st0
      fld qword [ftmp%]
      faddp st1,st0
      fistp dword [^xm%]          ;Calculate xm%
      
      fld qword [^Yd]             ;Calculate ym%
      fld st0
      fmulp st1,st0
      fld qword [^Md]
      fld st0
      fmulp st1,st0
      fsubrp st1,st0
      fsqrt
      fld qword [^Md]
      fdivp st1,st0
      
      fld qword [^R(3)]
      fmulp st1,st0
      fild dword [^rdl%]
      fmulp st1,st0
      fstp qword [ftmp%]
      
      fld qword [^ff]
      fld1
      fsubrp st1,st0
      fld qword [^ys]
      fmulp st1,st0
      fld qword [^ff]
      fld qword [^ye]
      fmulp st1,st0
      faddp st1,st0
      fild dword [^yscreen%]
      fmulp st1,st0
      fld qword [ftmp%]
      faddp st1,st0
      fistp dword [^ym%]          ;Calculate ym%
      
      call drw%
      
      inc dword [^c%]
      mov edx,[^c%]
      cmp edx,[^Lm%]
      jbe near loop
      
      ret
      ]
      endproc
      
      rem-----------------------------------------------------
      
      def proc_dump
      ff=1.0*c%/Lm%
      
      Zmm=Zs*(1.0-ff) + Ze*ff
      Fac=(1.0-A)/((Zmm+Dis)-A) : rem Reduction or scale factor
      
      xs=Xs*Fac : ys=Ys*Fac
      xe=Xe*Fac : ye=Ye*Fac
      
      rdl%=int(50.0-Zmm)
      
      xm%=int((xs*(1.0-ff) + xe*ff)*xscreen% + rdl%*R(1)*(sqr(Md*Md-Xd*Xd))/Md)
      ym%=int((ys*(1.0-ff) + ye*ff)*yscreen% + rdl%*R(3)*(sqr(Md*Md-Yd*Yd))/Md)
      endproc
      rem-----------------------------------------------------