3

Given the longitude and latitude coordinates of multiple points on a sphere, I can connect them using circular arc segments by using the Arc command in Asymptote:

path3 arc(triple c,real r,real theta1,real phi1,real theta2,real phi2,triple normal=O);

As an example, given (-10,0), (0,10), (10,0) and (-10,0) in degrees, I can use the following code to draw four black arcs on the sphere to make it look like I've connected the four points from start to finish (and cycled).

settings.render = 0;

import solids;
import three;
size(6cm);

real zenith  = pi/12.0;
real azimuth = pi/12.0;
currentprojection = perspective(cos(azimuth)*cos(zenith),
                                sin(azimuth)*cos(zenith),
                                sin(zenith));
real r = 1;
real ar = 1.2;

path3 myarc = Arc(c=O,normal=X, v1=-Z*r, v2=Z*r, n=12);
surface sphere = surface(myarc, angle1=0, angle2=360, c=O, axis=Z, n=12);

draw((0,0,  0)--(ar,0,0),   red+linewidth(1pt));
draw((0,0,  0)--(0,ar,0), green+linewidth(1pt));
draw((0,0,-ar)--(0,0,ar),  blue+linewidth(1pt));

label("$S$",(0,0,-ar),S);
draw(sphere, surfacepen=material(white+opacity(0.8),ambientpen=white));
label("$N$",(0,0,ar),N);

pair[] region = new pair[] {(-10,0),(0,10),(10,0),(0,-10),(-10,0)};
for(int i=1; i<region.length; ++i){
  pair one, two;
  if(region[i-1].y < region[i].y){
    one = (region[i]  .x, 90 - region[i]  .y);
    two = (region[i-1].x, 90 - region[i-1].y);
  }else if(region[i-1].y > region[i].y){
    one = (region[i-1].x, 90 - region[i-1].y);
    two = (region[i]  .x, 90 - region[i]  .y);
  }else if(region[i-1].x > region[i].x){
    one = (region[i]  .x, 90 - region[i]  .y);
    two = (region[i-1].x, 90 - region[i-1].y);
  }else{
    one = (region[i-1].x, 90 - region[i-1].y);
    two = (region[i]  .x, 90 - region[i]  .y);
  }
  path3 temp = arc(O,r,one.y,one.x,two.y,two.x);
  draw(temp,black+linewidth(.5pt));
}

Notes: the if - else if - else if - else block is used to draw the shorter one between the two possible arcs between each pair of point.

The code above produces this image without the 3D viewer being enabled.

The output of the code above

Now, I'd like to fill that region on the sphere with any color (for example, red with 50% transparency). How can I achieve this using Asymptote?

Improve this question
7
  • I'm aware that the sphere is completely covering the red axis, but I would like to think about it on my own for now. :)
    – Frenzy Li
    Jul 25, 2016 at 11:13
  • I have solved the issue in my first comment by drawing two hemispheres. I draw a hemisphere from -90 to 90 degrees, then I draw the red axis, after which I draw the second hemisphere from 90 to 270 degrees.
    – Frenzy Li
    Jul 25, 2016 at 11:49
  • You want to create a red surface covering the region you want to "fill." Since you're using settings.render=0, you need to draw the red surface after you've drawn the sphere. If you were using rasterized rendering, you'd instead make the red surface part of a sphere with a slightly larger radius (e.g., 1.001 instead of 1) so that the program knows which surface is closer to the viewer.
    – Charles Staats
    Jul 26, 2016 at 17:05
  • 1
    Additional note: since you're using arc and generating arcs of less than 90 degrees, I'm guessing each arc in the black patch consists of only one segment. You should be able to string these together into a cyclic path of length four using the & connector: path outline = arc1 & arc2 & arc3 & arc4 & cycle;. At this point, surface(patch(outline)) might give you a reasonable surface to fill the patch, at least with render=0.
    – Charles Staats
    Jul 26, 2016 at 17:08
  • @CharlesStaats Is it that: by drawing a cyclic path and patching it won't necessarily produce the spherical region, but an approximation to it that will land inside the sphere?
    – Frenzy Li
    Jul 26, 2016 at 17:12

1 Answer 1

Reset to default
2

Sorry I made a confusion in the first answer about the region (parallel circles and great circles), you can see the code at the end (the approximation is ok for small value of angle <30).

About the second solution, i.e. construct the region with a single patch, computing the derivative is easy (from the arc computation of Asymptote) but the obtained surface is far from the desired result (try to draw the oo3 surface in the code). I think that it is possible by decomposing the region into 4 or 8 parts but in any case it is crucial to compute the internal control points of the Bézier patch. It needs some knowledge about Bézier patch and surface approximation.

So I tried to complete the first solution: a parametric definition of the region. Since I do not want to define myself, I used the Arc definition of Asymptote. The strange parametrization (see the definition of f) is needed to avoid some numerical artefacts in the interpolation and you can recognize the Tchebychev points. It seems to be ok, for angle less or equal than 80. Please consider the code

import solids;
import three;
size(6cm);

real zenith  = pi/12.0;
real azimuth = pi/12.0;
currentprojection = perspective(cos(azimuth)*cos(zenith),
                                sin(azimuth)*cos(zenith),
                                sin(zenith));
real r = 1;
real ar = 1.2;

path3 myarc = Arc(c=O,normal=X, v1=-Z*r, v2=Z*r, n=12);
surface sphere = surface(myarc, angle1=0, angle2=360, c=O, axis=Z, n=12);

draw((0,0,  0)--(ar,0,0),   red+linewidth(1pt));
draw((0,0,  0)--(0,ar,0), green+linewidth(1pt));
draw((0,0,-ar)--(0,0,ar),  blue+linewidth(1pt));

label("$S$",(0,0,-ar),S);
draw(sphere, surfacepen=material(white+opacity(0.8),ambientpen=white));
label("$N$",(0,0,ar),N);
real angle=80;
path3[] aar;
pair[] region = new pair[] {(-angle,0),(0,angle),(angle,0),(0,-angle),(-angle,0)};
for(int i=1; i<region.length; ++i){
  pair one, two;
  if(region[i-1].y < region[i].y){
    one = (region[i]  .x, 90 - region[i]  .y);
    two = (region[i-1].x, 90 - region[i-1].y);
  }else if(region[i-1].y > region[i].y){
    one = (region[i-1].x, 90 - region[i-1].y);
    two = (region[i]  .x, 90 - region[i]  .y);
  }else if(region[i-1].x > region[i].x){
    one = (region[i]  .x, 90 - region[i]  .y);
    two = (region[i-1].x, 90 - region[i-1].y);
  }else{
    one = (region[i-1].x, 90 - region[i-1].y);
    two = (region[i]  .x, 90 - region[i]  .y);
  }
  path3 temp = Arc(O,r,one.y,one.x,two.y,two.x,32);
  draw(temp,black+linewidth(.5pt));
  aar.push(temp);
}
dot((cos(angle/180*pi),sin(angle/180*pi),0));
triple AA=(cos(angle/180*pi),sin(angle/180*pi),0);
dot((cos(angle/180*pi),sin(-angle/180*pi),0),blue);
triple BB=(cos(angle/180*pi),sin(-angle/180*pi),0);
dot((cos(angle/180*pi),0,-sin(angle/180*pi)),red);
triple CC=(cos(angle/180*pi),0,-sin(angle/180*pi));
dot((cos(angle/180*pi),0,sin(angle/180*pi)),green);
triple DD=(cos(angle/180*pi),0,sin(angle/180*pi));

real angr=angle/180*pi;
triple f(pair t)
{
  triple x1=point(aar[1],(1+cos((2*t.x+1)*pi))/2*length(aar[1]));
  triple x2=point(aar[3],(1+cos((2*t.x+1)*pi))/2*length(aar[3]));
  path3 temp_arc=Arc(O,x1,x2,32);
  return(point(temp_arc, (1+cos((2*t.y+1)*pi))/2*length(temp_arc)));
}
surface s=surface(f,(0,0),(1,1),32,32,Spline);
draw(shift(0.005*ar,0,0)*s,2bp+blue);

patch oo3=patch(AA{-dir(aar[1])}..{-dir(aar[1],0)}DD{dir(aar[0],0)}..{dir(aar[0])}BB{dir(aar[3],0)}..{dir(aar[3])}CC{-dir(aar[2])}..{-dir(aar[2],0)}cycle);

and the result (angle=80) with OpenGL renderer.

enter image description here

Previous anwser : With a few computations, it is possible to describe the region through a parametric surface and then to draw it.

In the spirit of Charles Staats indication, you can have an approximation of the surface using the patch constructor. By specifying the tangent at the vertices the approximation is ok for low angle. The best should be to compute the last four control points of the surface (see the definition of the sphere in three_surface.asy). To avoid some artefacts the surface is shifted and it works with render=0 and the OpenGL renderer.

import solids;
import three;
size(6cm);

real zenith  = pi/12.0;
real azimuth = pi/12.0;
currentprojection = perspective(cos(azimuth)*cos(zenith),
                               sin(azimuth)*cos(zenith),
                               sin(zenith));
real r = 1;
real ar = 1.2;

path3 myarc = Arc(c=O,normal=X, v1=-Z*r, v2=Z*r, n=12);
surface sphere = surface(myarc, angle1=0, angle2=360, c=O, axis=Z, n=12);

draw((0,0,  0)--(ar,0,0),   red+linewidth(1pt));
draw((0,0,  0)--(0,ar,0), green+linewidth(1pt));
draw((0,0,-ar)--(0,0,ar),  blue+linewidth(1pt));

label("$S$",(0,0,-ar),S);
draw(sphere, surfacepen=material(white+opacity(0.8),ambientpen=white));
label("$N$",(0,0,ar),N);

pair[] region = new pair[] {(-10,0),(0,10),(10,0),(0,-10),(-10,0)};
triple[] mypoints;
for(int i=1; i<region.length; ++i){
  pair one, two;
  if(region[i-1].y < region[i].y){
    one = (region[i]  .x, 90 - region[i]  .y);
    two = (region[i-1].x, 90 - region[i-1].y);
  }else if(region[i-1].y > region[i].y){
    one = (region[i-1].x, 90 - region[i-1].y);
    two = (region[i]  .x, 90 - region[i]  .y);
  }else if(region[i-1].x > region[i].x){
    one = (region[i]  .x, 90 - region[i]  .y);
    two = (region[i-1].x, 90 - region[i-1].y);
  }else{
    one = (region[i-1].x, 90 - region[i-1].y);
    two = (region[i]  .x, 90 - region[i]  .y);
  }
  path3 temp = arc(O,r,one.y,one.x,two.y,two.x);
  draw(temp,black+linewidth(.5pt));
}
real angle=10;
dot((cos(angle/180*pi),sin(angle/180*pi),0));
triple AA=(cos(angle/180*pi),sin(angle/180*pi),0);
dot((cos(angle/180*pi),sin(-angle/180*pi),0),blue);
triple BB=(cos(angle/180*pi),sin(-angle/180*pi),0);
dot((cos(angle/180*pi),0,-sin(angle/180*pi)),red);
triple CC=(cos(angle/180*pi),0,-sin(angle/180*pi));
dot((cos(angle/180*pi),0,sin(angle/180*pi)),green);
triple DD=(cos(angle/180*pi),0,sin(angle/180*pi));

patch oo=patch(AA{cross(AA,AA-CC)}..{cross(DD,DD-BB)}DD{cross(DD,DD-AA)}.. {cross(BB,BB-CC)}BB{cross(BB,BB-DD)}..{cross(CC,CC-AA)}CC{cross(CC,AA-DD)}..{cross(AA,AA-DD)}cycle);
draw(shift(.005*ar,0,0)*surface(oo),blue);

And the result enter image description here

For angle=45, the difference is apparent (with the OpenGL renderer)...

Improve this answer
3
  • This actually works pretty well for small angles, but as you say, angle=45 or angle=90, the patch doesn't have the correct tangents. I think if the tangent leaving AA would be in the same plane with the great circle determined by AA and BB, the difference will be gone.
    – Frenzy Li
    Aug 4, 2016 at 13:25
  • Not sure that the tangents are wrong. A Bezier patch has 16 control points, here 12 are defined and the last 4 are determined by Asymptote. You can observe a small difference between the two definitions real a=4/3*(sqrt(2)-1);patch oo1=patch(X{Y}..{-X}Y{Z}..{-Y}Z..Z{X}..{-Z}cycle, new triple[] {(1,a,a),(a,1,a),(a^2,a,1), (a,a^2,1)}); patch oo2=patch(X{Y}..{-X}Y{Z}..{-Y}Z..Z{X}..{-Z}cycle); oo1 is octant defined in three_surface.asy to construct the sphere.
    – O.G.
    Aug 4, 2016 at 15:28
  • @Frenzy DT : You're right, the tangents should be modified.
    – O.G.
    Aug 4, 2016 at 16:54

You must log in to answer this question.

Not the answer you're looking for? Browse other questions tagged .