Subversion Repositories wimsdev

Rev

Rev 18320 | Rev 18326 | Go to most recent revision | Show entire file | Ignore whitespace | Details | Blame | Last modification | View Log | RSS feed

Rev 18320 Rev 18325
Line 648... Line 648...
648
    if (tikz_file)
648
    if (tikz_file)
649
      fprintf(tikz_file, "(%i,%i)--(%i,%i)",
649
      fprintf(tikz_file, "(%i,%i)--(%i,%i)",
650
       pm->p[i],flip(pm->p[i+1]),pm->p[i+2],flip(pm->p[i+3]));
650
       pm->p[i],flip(pm->p[i+1]),pm->p[i+2],flip(pm->p[i+3]));
651
  }
651
  }
652
  if(tikz_file) fprintf(tikz_file, ";\n");
652
  if(tikz_file) fprintf(tikz_file, ";\n");
-
 
653
}
-
 
654
/* for geodesics, return data for obj_hypgeods + a point on the arc
-
 
655
  to be used for filling hyperbolic triangle*/
-
 
656
int hypgeodaux(double *q, double* res)
-
 
657
{
-
 
658
  double r,cx,cy,a1,a2,a3,tmp,
-
 
659
    nx = -q[0]*q[2]*q[2]+(q[0]*q[0]+q[1]*q[1]+1)*q[2]-q[0]*q[3]*q[3]-q[0],
-
 
660
    ny = -q[1]*q[2]*q[2]-q[1]*q[3]*q[3]+(q[0]*q[0]+q[1]*q[1]+1)*q[3]-q[1],
-
 
661
    dy = -2*q[1]*q[2]+2*q[0]*q[3];
-
 
662
  if (dy*dy*1e4 < nx*nx+ny*ny){
-
 
663
    res[5]=(q[0]+q[2])/2;
-
 
664
    res[6]=(q[1]+q[3])/2;
-
 
665
    return 0;}
-
 
666
  cx = ny/dy; cy=-nx/dy;
-
 
667
  a1 = atan2(q[1]-cy, q[0]-cx);
-
 
668
  a2 = atan2(q[3]-cy, q[2]-cx);
-
 
669
  a3 = (a1+a2)/2;
-
 
670
  if (fabs(a2-a1)>M_PI){if(a1<a2) a1+=2*M_PI; else a2+=2*M_PI;};
-
 
671
  if(a1>a2) {tmp=a1; a1=a2; a2=tmp;}
-
 
672
  r = sqrt(cx*cx+cy*cy-1);
-
 
673
  res[0]=cx;
-
 
674
  res[1]=cy;
-
 
675
  res[2]=r;
-
 
676
  res[3]=a1;
-
 
677
  res[4]=a2;
-
 
678
  res[5]=cx+r*cos(a3);
-
 
679
  res[6]=cy+r*sin(a3);
-
 
680
  return 1;
653
}
681
}
654
/* hyperbolic geodesics from x1,y1 to x2,y2, x3,y3 to x4,y4 in Poincaré disk */
682
/* hyperbolic geodesics from x1,y1 to x2,y2, x3,y3 to x4,y4 in Poincaré disk */
655
 
683
 
656
void obj_hypgeods(objparm *pm)
684
void obj_hypgeods(objparm *pm)
657
{
685
{
658
 int i, *qq = pm->p;
686
  int i, *qq=pm->p;
-
 
687
  double rx,ry,res[7];
659
 for (i = 0; i < pm->pcnt; i+=4){
688
  for (i = 0; i < pm->pcnt; i+=4)
660
   double r, *q=pm->pd + i,
689
    if (hypgeodaux(pm->pd+i,res)){
-
 
690
      scale2(res[2],res[2],&rx,&ry);
-
 
691
      scale(res,qq,1);
661
   nx = -q[0]*q[2]*q[2]+(q[0]*q[0]+q[1]*q[1]+1)*q[2]-q[0]*q[3]*q[3]-q[0],
692
      myGdImageArc(image,qq[0],qq[1],rx,ry,res[3]/DEG,res[4]/DEG,pm->color[0]);
-
 
693
    }
-
 
694
    else {
-
 
695
      scale(pm->pd+i,qq,2);
662
   ny = -q[1]*q[2]*q[2]-q[1]*q[3]*q[3]+(q[0]*q[0]+q[1]*q[1]+1)*q[3]-q[1],
696
      gdImageLine(image,qq[0],qq[1],qq[2],qq[3],pm->color[0]);
-
 
697
    }
-
 
698
}
-
 
699
/* hyperbolique triangle, can be filled */
-
 
700
void obj_hyptriangle(objparm *pm)
-
 
701
{
663
   dy = -2*q[1]*q[2]+2*q[0]*q[3];
702
  double rx,ry,data[8],res[7];
-
 
703
  int i, *qq=pm->p;
-
 
704
  for(i=0; i<8; i++) data[i]=pm->pd[i%6];
-
 
705
  for(i=0; i<3; i++)
664
   if (dy*dy * 1e8 < nx*nx+ny*ny){
706
    if(hypgeodaux(data+2*i,res)) {
-
 
707
      scale2(res[2],res[2],&rx,&ry);
665
     scale(q,qq,2);
708
      scale(res,qq,1);
-
 
709
      myGdImageArc(image,qq[0],qq[1],rx,ry,res[3]/DEG,res[4]/DEG,pm->color[0]);
-
 
710
    }
-
 
711
    else {
-
 
712
      scale(data+2*i,qq,2);
666
     gdImageLine(image,qq[0],qq[1],qq[2],qq[3],pm->color[0]);
713
      gdImageLine(image,qq[0],qq[1],qq[2],qq[3],pm->color[0]);
-
 
714
    }
-
 
715
  if(pm->fill){
-
 
716
    data[0]=res[5];
-
 
717
    data[1]=res[6];
-
 
718
    hypgeodaux(data,res);
-
 
719
    scale(res+5,qq,1);
-
 
720
    patchgdImageFill(image,qq[0],qq[1],pm->color[0]);
-
 
721
  }
-
 
722
}
-
 
723
/* hyperbolic geodesic through two points */
-
 
724
void obj_hypline(objparm *pm)
-
 
725
{
-
 
726
  double rx,ry,res[7];
-
 
727
  if (hypgeodaux(pm->pd,res)){
-
 
728
      double alpha=atan(1/res[2]), beta = atan2(res[0], res[1]);
-
 
729
      scale(pm->pd, pm->p, 1);
-
 
730
      scale2(res[2],res[2],&rx,&ry);
-
 
731
      myGdImageArc(image, pm->p[0], pm->p[1], rx, ry, 180+(beta-alpha)/DEG, 180+(beta+alpha)/DEG,pm->color[0]);
-
 
732
  }
-
 
733
  else {
-
 
734
    double gamma = atan2(pm->pd[0],pm->pd[1]);
-
 
735
    res[0]=cos(gamma);
-
 
736
    res[1]=sin(gamma);
-
 
737
    res[2]=-cos(gamma);
-
 
738
    res[3]=-sin(gamma);
-
 
739
    scale(res,pm->p,2);
-
 
740
    gdImageLine(image,pm->p[0],pm->p[1],pm->p[2],pm->p[3],pm->color[0]);
-
 
741
  }
-
 
742
}
-
 
743
 
-
 
744
void obj_hyplines(objparm *pm)
-
 
745
{
-
 
746
  double rx,ry,res[7];
-
 
747
  int i;
-
 
748
  for (i=0; i < pm->pcnt; i+=4)
-
 
749
    if (hypgeodaux(pm->pd+i,res)) {
-
 
750
      double alpha=atan(1/res[2]), beta = atan2(res[0], res[1]);
-
 
751
      scale(pm->pd+i, pm->p, 1);
-
 
752
      scale2(res[2],res[2],&rx,&ry);
-
 
753
      myGdImageArc(image, pm->p[0], pm->p[1], rx, ry, 180+(beta-alpha)/DEG, 180+(beta+alpha)/DEG,pm->color[0]);
667
   }
754
    }
668
    else {
755
    else {
669
      double rx,ry;
-
 
670
      double dc[2];int ic[2];
-
 
671
      double
-
 
672
        cx = ny/dy, cy=-nx/dy,
-
 
673
        a1 = atan2(q[1]-cy, q[0]-cx)/DEG,
756
      double gamma = atan2(pm->pd[i],pm->pd[i+1]);
674
        a2 = atan2(q[3]-cy, q[2]-cx)/DEG;
757
      res[0]=cos(gamma);
675
      if (fabs(a2-a1)>180){if(a1<a2) a1+=360; else a2+=360;};
-
 
676
      r = sqrt(cx*cx+cy*cy-1);
758
      res[1]=sin(gamma);
677
      scale2(r,r,&rx,&ry);
759
      res[2]=-cos(gamma);
678
      dc[0]=cx;dc[1]=cy;
760
      res[3]=-sin(gamma);
679
      scale(dc,ic,1);
761
      scale(res,pm->p,2);
680
      if(a1<a2) myGdImageArc(image,ic[0],ic[1],rx,ry,a1,a2,pm->color[0]);
-
 
681
      else myGdImageArc(image,ic[0],ic[1],rx,ry,a2,a1,pm->color[0]);
762
      gdImageLine(image,pm->p[0],pm->p[1],pm->p[2],pm->p[3],pm->color[0]);
682
    }
763
    }
683
  }
-
 
684
}
764
}
685
/* segments */
765
/* segments */
686
void obj_dlines(objparm *pm)
766
void obj_dlines(objparm *pm)
687
{
767
{
688
 int i, n;
768
 int i, n;
Line 746... Line 826...
746
{
826
{
747
  int i,nb;
827
  int i,nb;
748
  double dt, r=rx>ry?rx:ry;
828
  double dt, r=rx>ry?rx:ry;
749
  t2-=t1;
829
  t2-=t1;
750
  t2-=360*floor(t2/360);
830
  t2-=360*floor(t2/360);
751
  if(t2<1) t2=360;
831
  if(t2==0) t2=360;
752
  t1*=DEG;t2*=DEG;
832
  t1*=DEG;t2*=DEG;
753
  nb = r*t2/3;
833
  nb = r*t2/3;
754
  dt = t2/nb;
834
  dt = t2/nb;
755
  for (i=0; i<nb; ++i,t1+=dt)
835
  for (i=0; i<nb; ++i,t1+=dt)
756
    gdImageLine(im,cx+rx*cos(t1),cy+ry*sin(t1),
836
    gdImageLine(im,cx+rx*cos(t1),cy+ry*sin(t1),