// Grafica un diagrama en hojas de un paraboloide
// parábola: y= x^2/4f

// datos iniciales
// k = cantidad de hojas, f = distancia focal, 
// largo= nº de puntos para graficar la longitud de la hoja,
// tm = valor nominal de la longitud de una hoja
// tam = tamaño en pixeles de la gráfica
//rhom = valor nominal del radio máximo del paraboloide
k=16; f=0.35; largo =501; lar= ceil(6.3*largo/k); tam= 500; rhom=2.35;

//coordenadas
rho = linspace(0.0001,1, largo); rho= rho.*rho; rho= rhom*rho
L= 0.5*rho.*sqrt(1+0.25/f^2*rho.*rho) + f*log(0.5/f*rho+sqrt(1+0.25/f^2*rho.*rho))
alfa = %pi/k*(rho./L)
x=L.*cos(alfa)
y=L.*sin(alfa)
ang=alfa(largo)*linspace(1,0, lar)
xs= L(largo)*cos(ang)
ys= L(largo)*sin(ang)
x=[x,xs,xs($:-1:1),x($:-1:1)];
y=[y,ys,-ys($:-1:1),-y($:-1:1)];
//u=[x,xs]; v=[y,ys]
//plot2d(u,v)
coor=[x;y];
long=L(largo)

//preparación gráfica
a=gca();
a.isoview="on";
a.parent.figure_size= [tam,tam];
//xsetech(wrect=[0,0,1.0,1.0],arect=1/60*[1,1,1,1])
//square(-long, -long, long, long)

//rotación y graficación
for i=2*%pi*(1:k)/k
xy=rotate(coor,i)
plot2d(xy(1,:),xy(2,:), axesflag=0)
end

//graficación del cuadrado (diseños optimizado, 
//con lámina superpuesta al diagrama de hojas)
u=long*[-1, 0, 1, 0; 0, 1, 0, -1]
v=long*[0, 1, 0, -1; 1, 0, -1, 0]
xsegs(u,v)

nomar='fig3' //nombre del archivo gráfico de salida
//exportar a PDF
//scf(0); xs2pdf(0,nomar); xs2pdf(gcf(),nomar);
//exportar a SVG
//scf(0); xs2svg(0,nomar); xs2svg(gcf(),nomar);