% hyperboloïde à une nappe
% modèle d'une tour de refroidissement
% x=a*sqrt(1+u^2)*cos(v)
% y=a*sqrt(1+u^2)*sin(v)
% z=c*u+H
% -H/c < u < 0

/HyperboloidDict  100 dict def
HyperboloidDict begin
/MaillageHyperboloid{
 0 incrementHAUTEUR H incrementHAUTEUR 0.8 mul sub {
 /z1 ED
 /u1 {z1 H sub c div} bind def
 /r1 {1 u1 dup mul add sqrt a mul} bind def
0 dT 360 {
    /v ED
% M1
    /x1 r1 v cos mul def
    /y1 r1 v sin mul def
    /Xpoint x1 def
    /Ypoint y1 def
    /Zpoint z1 def
PointsApresTransformations
    /Z1 ED
    /Y1 ED
    /X1 ED
% M2
    /x2 r1 v dT add cos mul def
    /y2 r1 v dT add sin mul def
    /Xpoint x2 def
    /Ypoint y2 def
    /Zpoint z1 def
PointsApresTransformations
    /Z2 ED
    /Y2 ED
    /X2 ED
% M3
    /z3 z1 incrementHAUTEUR add def
    /u3 {z3 H sub c div} bind def
    /r3 {1 u3 dup mul add sqrt a mul} bind def
    /x3 r3 v dT add cos mul def
    /y3 r3 v dT add sin mul def
    /Xpoint x3 def
    /Ypoint y3 def
    /Zpoint z3 def
PointsApresTransformations
    /Z3 ED
    /Y3 ED
    /X3 ED
% M4
    /x4 r3 v cos mul def
    /y4 r3 v sin mul def
    /Xpoint x4 def
    /Ypoint y4 def
    /Zpoint z3 def
PointsApresTransformations
    /Z4 ED
    /Y4 ED
    /X4 ED
%% centre de la face
   /xF {X1 X2 add X3 add X4 add 4 div} bind def
   /yF {Y1 Y2 Y3 Y4 add add add 4 div} bind def
   /zF {Z1 Z2 Z3 Z4 add add add 4 div} bind def
%% vecteur : centre de la face -> point de vue
    /xV {XpointVue xF sub } bind def
    /yV {YpointVue yF sub } bind def
    /zV {ZpointVue zF sub } bind def
%% distance
%    /D {xV dup mul yV dup mul add zV dup mul add sqrt} bind def
%% normale
% les vecteurs de la face
 /V1x X2 X1 sub def
 /V1y Y2 Y1 sub def
 /V1z Z2 Z1 sub def
 /V2x X3 X2 sub def
 /V2y Y3 Y2 sub def
 /V2z Z3 Z2 sub def
% le produit vectoriel
    /xN V1y V2z mul V1z V2y mul sub def
    /yN V1x V2z mul V1z V2x mul sub neg def
    /zN V1x V2y mul V1y V2x mul sub def
% le produit scalaire des 2 vecteurs
/PSfacette xV xN mul yV yN mul add zV zN mul add def
% vecteur centre de la facette vers la source de lumière
    /VxLight xLight xF sub def
    /VyLight yLight yF sub def
    /VzLight zLight zF sub def
% norme
/normeL  VxLight dup mul
         VyLight dup mul
         VzLight dup mul
         add add sqrt def
% produit scalaire
/PSLight  VxLight xN mul
          VyLight yN mul add
          VzLight zN mul add def
/normeN xN dup mul yN dup mul zN dup mul add add sqrt def
/CosCouleur PSLight normeL normeN mul div ChangeSigne def
condition {
/facette{
newpath
X1 Y1 Z1 3dto2d moveto
X2 Y2 Z2 3dto2d lineto
X3 Y3 Z3 3dto2d lineto
X4 Y4 Z4 3dto2d lineto
closepath
} def
facette
nocolor
toDraw} if
 }for
    }for
} def
%
/Bases {
% centre de la base inférieure après transformations
    /Xpoint 0 def
    /Ypoint 0 def
    /Zpoint 0 def
  PointsApresTransformations
    /zCinf ED
    /yCinf ED
    /xCinf ED 
% centre de la base supérieure avant transformations
    /Xpoint 0 def
    /Ypoint 0 def
    /Zpoint H def
  PointsApresTransformations
    /zCsup ED
    /yCsup ED
    /xCsup ED
% Normale extérieure à la base inférieure
    /nXBaseInf xCinf xCsup sub def
    /nYBaseInf yCinf yCsup sub def
    /nZBaseInf zCinf zCsup sub def 
% Normale extérieure à la base supérieure
    /nXBaseSup xCsup xCinf sub def
    /nYBaseSup yCsup yCinf sub def
    /nZBaseSup zCsup zCinf sub def 
% rayon vers point de vue
    /RXvueSup XpointVue xCsup sub def
    /RYvueSup YpointVue yCsup sub def
    /RZvueSup ZpointVue zCsup sub def
    /RXvueInf XpointVue xCinf sub def
    /RYvueInf YpointVue yCinf sub def
    /RZvueInf ZpointVue zCinf sub def
% Visibilité de la base inférieure
    /PSbaseInf nXBaseInf RXvueInf mul
               nYBaseInf RYvueInf mul add
               nZBaseInf RZvueInf mul add def
% Visibilité de la base supérieure
    /PSbaseSup nXBaseSup RXvueSup mul
               nYBaseSup RYvueSup mul add
               nZBaseSup RZvueSup mul add def }
def

/FaceInf {
/TableauxPoints [
0 1 359 {% on décrit le cercle
    /theta exch def [
    /Xpoint R theta cos mul def
    /Ypoint R theta sin mul def
    /Zpoint 0 def
  PointsApresTransformations
   3dto2d ]
    } for
    ] def
 TableauxPoints 0 get aload pop moveto
0 1 359 {
    /compteur exch def
    TableauxPoints compteur get aload pop
    lineto } for }
def


/FaceSup {
/TableauxPoints [
0 1 359 {% on décrit le cercle
    /theta exch def [
    /Xpoint r theta cos mul def
    /Ypoint r theta sin mul def
    /Zpoint H def
  PointsApresTransformations
   3dto2d ]
    } for
    ] def
 TableauxPoints 0 get aload pop moveto
0 1 359 {
    /compteur exch def
    TableauxPoints compteur get aload pop
    lineto } for }
 def
end