外部変形は データのやり取りをテキストファイルで行うので プログラム言語は 自由に選ぶことができます。図形は機能的かつシンプルなため、数多くのユーザーに受け入れられています。
gccでソリッド図形のＨ形鋼断面を描く
:gccでソリッド図形のＨ形鋼断面を描く
@echo off
if not exist %~dpn0.exe (
for /f "delims=:" %%n in ('findstr /n "^#!" %0') do (
more +%%n %0 > %~dpn0.c
gcc -Os %~dpn0.c -o %~dpn0.exe -s
)
)
set s=%*
if defined s (
echo ^%*> %~dpn0.txt
) else (
if exist %~dpn0.txt (
for /f "tokens=*" %%a in (%~dpn0.txt) do set s=%%a
)
)
if not defined s set s=400 200 8 13 16
%~dpn0
goto:eof
REM #jww
REM #1%d 位置を指示して下さい
REM #99#
REM #c H B T1 T2 R =
REM #e
#!この次の行からプログラムを書いてください
#include "jw.h"
void slplot(FILE *f, double *p1, double *p2, double *p3, double *p4)
{
fprintf(f, "sl %.15g %.15g %.15g %.15g %.15g %.15g %.15g %.15g\n", p1[0], p1[1], p2[0], p2[1], p3[0], p3[1], p4[0], p4[1]);
}
void jish(FILE *f, double *hp, double h, double b, double tw, double tf, double r, double d, int pt)
{
double p0[2], pa[2], pb[2], pc[2], x0, y0, x1, y1, drad, x, y, a;
double p[5][2];
double pp[20][2];
int u = 5;
int i, j, k;
int red = 0;
int green = 128;
int blue = 128;
fprintf(f, "lc10 %d", red + 256 * green + 256 * 256 * blue);
fprintf(f, "pl\n");
offset(b, h, pt, 5, p0);
x0 = p0[0];
y0 = p0[1];
x1 = hp[0];
y1 = hp[1];
drad = rad(d);
p[0][0] = tw / 2; p[0][1] = h / 2 - tf - r*0;
p[1][0] = tw / 2 + r; p[1][1] = h / 2 - tf - r;
p[2][0] = tw / 2 + r; p[2][1] = h / 2 - tf;
p[3][0] = b / 2; p[3][1] = h / 2 - tf;
p[4][0] = b / 2; p[4][1] = h / 2;
for (j = 0; j <= 3; j++) {
a = 90 * ( (1 - j) % 4);
for (i = 0; i < u; i++) {
if (j % 2 == 0) { k = i; } else { k = u - i - 1; }
x = p[k][0];
y = p[k][1];
if (j % 3 > 0) x = -x;
if (j > 1) y = -y;
pa[0] = x + x0;
pa[1] = y + y0;
rot(pa, drad, pb);
moveto(pb, x1, y1, pc);
if (k == 1 && r > 0){
fprintf(f, "sc %.15g %.15g %.15g %.15g %.15g %.15g %.15g -1\n", pc[0], pc[1], r, 1.0, drad, rad(a), PI / 2);
}
pp[i + j * u][0] = pc[0];
pp[i + j * u][1] = pc[1];
}
}
slplot(f, pp[0], pp[9], pp[10], pp[19]);
slplot(f, pp[3], pp[4], pp[5], pp[6]);
slplot(f, pp[13], pp[14], pp[15], pp[16]);
fprintf(f, "#\n");
}
int main(void)
{
FILE *f;
char S_[256], *F[20];
double hk;
int i, n, hpn=0;
char nn[4];
double hp[101][2];
int pt = 5;
double h, b, t1, t2, r;
split(getenv("s"), F);
h = atof(F[0]);
b = atof(F[1]);
t1 = atof(F[2]);
t2 = atof(F[3]);
r = atof(F[4]);
if ( (f = fopen("jwc_temp.txt", "r")) != NULL) {
while (fgets(S_, 256, f) != NULL) {
split(chomp(S_), F);
if (strncmp(S_, "hk", 2) == 0) hk = atof(F[1]);
if (strncmp(S_, "hp", 2) == 0) {
n = strlen(F[0]);
substr(F[0], 2, n-1, nn);
chomp(nn);
if (n > 0 && nn[n-1] == '-')
{
substr(F[0], 2, n-2, nn);
}
i = atoi(nn);
hp[i][0] = atof(F[1]);
hp[i][1] = atof(F[2]);
hpn++;
}
}
} else {
exit(MISSING_JWC_TEMP_TXT);
}
fclose(f);
f = fopen("jwc_temp.txt", "w");
for (i = 1; i <= hpn; i++) {
jish(f, hp[i], h, b, t1, t2, r, hk, pt);
}
fclose(f);
return 0;
}
○jw.h
/* jw.h : JW_CAD 外部変形 例題用ヘッダファイル
動作確認
MinGW + gcc / g++ 4.5.0
Borland C++ 5.5.1 for Win32
Visual C++ 2005 Express Edition
*/
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>
#define PI 3.1415926535897932384626433832795
#define MISSING_JWC_TEMP_TXT printf("jwc_temp.txt が 見当たりません\n")
#define NO_OPERATION printf("処理できませんでした\n")
double arc_l(double r, double p1, double p2, double w);
double __ee__[151], __hh__[151];
int split(char *S_, char **F)
{ /* char S_[256], *F[20]; */
int i, NF;
for (i = 0; i < 20; i++) {
if ( (F[i] = strtok(S_, " ")) == NULL) { NF = i; break; }
S_ = NULL;
}
return NF;
}
char *substr(char *src, int off, int len, char *buf)
{
return strncpy(buf, src+off, len);
}
char *chomp(char *str)
{
int l = strlen(str);
if ( l > 0 && str[l-1] == '\n' )
{
str[l-1] = '\0';
}
return str;
}
int __sgn(double x)
{
if (x == 0) { return 0;
} else if (x > 0) { return 1;
} else { return -1;
}
}
double __abs(double x)
{
return (x > 0) ? x : -x;
}
int __int(double a)
{
int b;
if (a < 0) { b = ceil(a); } else { b = floor(a); }
return b;
}
int angle_eql(double d1, double d2, double eps)
{
double dd = __abs(d1 - d2);
while (dd > PI - eps) dd -= PI;
return __abs(dd) < eps ? 0 : 1;
}
int length_eql(double l1, double l2, double eps)
{
return __abs(l1 - l2) < eps ? 0 : 1;
}
double deg(double x)
{
return x * 180.0 / PI;
}
double rad(double x)
{
return x * PI / 180.0;
}
double hypot(double x, double y)
{
return sqrt(x * x + y * y);
}
int htoi(char *x)
{
int y;
if (strcmp(x, "a") == 0) { y = 10;
} else if (strcmp(x, "b") == 0) { y = 11;
} else if (strcmp(x, "c") == 0) { y = 12;
} else if (strcmp(x, "d") == 0) { y = 13;
} else if (strcmp(x, "e") == 0) { y = 14;
} else if (strcmp(x, "f") == 0) { y = 15;
} else {
y = atoi(x);
}
return y;
}
void ln_reset(double *ln, double *p)
{
double x1, y1, x2, y2, xs, ys;
x1 = ln[0];
y1 = ln[1];
x2 = ln[2];
y2 = ln[3];
/* 指示線の角度を -PI/2 < d <= +PI/2 とする */
if (x1 > x2) {
xs = x1; x1 = x2; x2 = xs;
ys = y1; y1 = y2; y2 = ys;
}
if (x1 == x2 && y1 > y2) {
ys = y1; y1 = y2; y2 = ys;
}
p[0] = x1;
p[1] = y1;
p[2] = x2;
p[3] = y2;
}
void ci_reset(double *ci, double *p)
{
double r, p1, p2, w, d;
r = ci[2];
p1 = ci[3] ? ci[3] : 0;
p2 = ci[4] ? ci[4] : 360;
w = ci[5] ? ci[5] : 1;
d = ci[6] ? ci[6] : 0;
if (w > 1) {
w = 1.0 / w;
r /= w;
p1 -= 90.0;
p2 -= 90.0;
d += 90.0;
}
if (p1 == p2) { p1 = 0; p2 = 0; }
if (__abs(p1 - p2) < 1e-12) { p2 += 360; }
/* if (__abs(p1 - p2) < 1e-12) { p1 = 0; p2 = 360; d = 0; } */
p[0] = ci[0];
p[1] = ci[1];
p[2] = r;
p[3] = p1;
p[4] = p2;
p[5] = w;
p[6] = d;
}
void polar(double r, double d, double w, double *p) {
p[0] = r * cos(d);
p[1] = r * sin(d) * w;
}
void rot(double *x, double d, double *p) {
double x1, y1, co, si;
x1 = x[0];
y1 = x[1];
co = cos(d);
si = sin(d);
p[0] = x1 * co - y1 * si;
p[1] = x1 * si + y1 * co;
}
void offset(double w, double h, int pt, int opt, double *p) {
p[0] = w * ( (9 - pt) % 3 - (9 - opt) % 3) / 2;
p[1] = h * (__int( (9 - pt) / 3) - __int( (9 - opt) / 3)) / 2;
}
void moveto(double *pt, double x, double y, double *p)
{
p[0] = pt[0] + x;
p[1] = pt[1] + y;
}
void polarto(double *pt, double r, double d, double *p)
{
p[0] = pt[0] + r * cos(d);
p[1] = pt[1] + r * sin(d);
}
double ptdist(double *p1, double *p2)
{
return sqrt(pow(p2[0] - p1[0], 2.0) + pow(p2[1] - p1[1], 2.0));
}
double ptslope(double *p1, double *p2)
{ /* ２点の角度を返す */
return atan2(p2[1] - p1[1], p2[0] - p1[0]);
}
void ptspan(double *p1, double *p2, double d, double *p)
{ /* d : 角度rad */
double co, si, x1, y1, x2, y2;
si = sin(d);
co = cos(d);
x1 = p1[0];
y1 = p1[1];
x2 = p2[0];
y2 = p2[1];
p[0] = (x2 - x1) * co + (y2 - y1) * si;
p[1] = - (x2 - x1) * si + (y2 - y1) * co;
}
void ptcenter(double *p1, double *p2, double x, double *p)
{
polarto(p1, ptdist(p1, p2) * x, ptslope(p1, p2), p);
}
int ptsgn(double *x, double *y, double d)
{ /* ２点を結ぶ線の傾斜の判別 ／ : 1, ＼ : -1 */
double p[2];
ptspan(x, y, 0, p);
rot(p, rad(-d), p);
return __sgn(p[0] * p[1]);
}
void lnpoint(double *ln, int pm, double *p)
{
if (pm == -1) {
p[0] = ln[0];
p[1] = ln[1];
} else if (pm == 1) {
p[0] = ln[2];
p[1] = ln[3];
} else {
p[0] = (ln[0] + ln[2]) / 2;
p[1] = (ln[1] + ln[3]) / 2;
}
}
double lnlength(double *ln)
{ /* 線長を返す */
return sqrt(pow(ln[2] - ln[0], 2.0) + pow(ln[3] - ln[1], 2.0));
}
double lnslope(double *ln)
{ /* 線の角度を返す */
return atan2(ln[3] - ln[1], ln[2] - ln[0]);
}
void lnspan(double *ln, double d, double *p)
{ /* 線の dx, dy を返す */
/* d : 角度rad */
double co, si, x1, y1, x2, y2;
si = sin(d);
co = cos(d);
x1 = ln[0];
y1 = ln[1];
x2 = ln[2];
y2 = ln[3];
p[0] = (x2 - x1) * co + (y2 - y1) * si;
p[1] = - (x2 - x1) * si + (y2 - y1) * co;
}
void lncenter(double *s, double x, double *p)
{
lnpoint(s, -1, p);
polarto(p, lnlength(s) * x, lnslope(s), p);
}
void lnmove(double *ln, double b, double pmy, double *l)
{ /* 線 ln を直交距離 b だけ pmy の方向に平行移動する */
double x1, y1, x2, y2, d, co, si;
x1 = ln[0];
y1 = ln[1];
x2 = ln[2];
y2 = ln[3];
d = lnslope(ln);
co = b * cos(d) * pmy;
si = b * sin(d) * pmy;
l[0] = x1 - si;
l[1] = y1 + co;
l[2] = x2 - si;
l[3] = y2 + co;
}
double ln_hpdist(double *ln, double *hp)
{ /* 線と点の最短距離を返す */
double x[4];
x[0] = hp[0];
x[1] = hp[1];
x[2] = ln[0];
x[3] = ln[1];
return lnlength(x) * sin(lnslope(ln) - lnslope(x));
}
void ln_hpsnap(double *ln, double *hp, double *p)
{ /* 線の吸着点を返す */
double t, d;
t = ln_hpdist(ln, hp);
d = lnslope(ln);
p[0] = hp[0] + t * sin(d);
p[1] = hp[1] - t * cos(d);
}
int ln_hpends(double *ln, double *hp)
{
double p1[2], p2[2], l1, l2;
lnpoint(ln,-1, p1);
lnpoint(ln, 1, p2);
l1 = ptdist(hp, p1);
l2 = ptdist(hp, p2);
return l1 < l2 ? -1 : 1;
}
int ln_hpsens(double *ln, double *hp)
{
double p1[2], p2[2], p4[2], ll, l1, l2;
lnpoint(ln,-1, p1);
lnpoint(ln, 1, p2);
ln_hpsnap(ln, hp, p4);
ll = lnlength(ln);
l1 = ptdist(p1, p4);
l2 = ptdist(p2, p4);
return (l1 + l2 == ll && (l1 == ll || l2 == ll)) ? 0 : ( (l1 <= ll && l2 <= ll) ? -1 : 1); /* -1 : 内, 1 : 外 */
}
int ln_hpmens(double *s, double *p)
{
return __sgn(ln_hpdist(s, p));
}
void cipoint(double *ci, int pm, double *p)
{
double x, y, r, p1, p2, w, d;
double co, si, xr, yr, q;
x = ci[0];
y = ci[1];
r = ci[2];
p1 = ci[3] ? ci[3] : 0;
p2 = ci[4] ? ci[4] : 360;
w = ci[5] ? ci[5] : 1;
d = ci[6] ? ci[6] : 0;
d *= PI / 180;
co = cos(d) * r;
si = sin(d) * r;
if (pm == -1) { q = p1 * PI / 180; } else { q = p2 * PI / 180; }
xr = cos(q);
yr = sin(q) * w;
p[0] = x + xr * co - yr * si;
p[1] = y + xr * si + yr * co;
p[2] = atan2(yr / w, xr * w) + d;
}
int cicenter(double c, double x, double *p) {
double nl, r, p1, p2, w, d;
double pc[2], cc, ep, pp, d1, d2, dd, __pl, __ml, l, eps, pr[2], ps[2], el;
int wflg;
pc[0] = c[0];
pc[1] = c[1];
r = c[2];
p1 = c[3] ? c[3] : 0;
p2 = c[4] ? c[4] : 360;
w = c[5] ? c[5] : 1;
d = c[6] ? c[6] : 0;
wflg = 0;
if (w > 1.0) {
w = 1.0 / w;
r /= w;
p1 -= 90.0;
p2 -= 90.0;
d += 90.0;
wflg = 1;
}
nl = arc_l(r, p1, p2, w) * x;
if (x == 0) {
cipoint(c, -1, p);
return 0;
}
cc = 10;
ep = r / pow(10, cc); /* 誤差 */
if (p1 == p2) p2 += 360;
/* if (__abs(p1 - p2) < 1e-12) { p1 = 0; p2 = 360; } */
pp = p2 - p1;
if (pp <= 0) pp += 360;
d1 = p1; d2 = p1; dd = pp * x;
__pl = nl + ep;
__ml = nl - ep;
if (w != 1) {
l = arc_l(1, 0, dd, w);
dd *= (rad(dd) / l);
} else {
l = nl;
}
eps = dd;
d2 += eps;
/* はさみうち法（ここから）*/
if (w != 1) {
l = arc_l(r, d1, d2, w);
while (l < __ml || __pl < l) {
eps /= 2.0;
if (eps < ep / r) break;
if (l < nl) { d2 += eps; } else { d2 -= eps; }
l = arc_l(r, d1, d2, w);
}
}
/* はさみうち法（ここまで）*/
if (l != 0) {
polar(r, rad(d2), w, pr);
rot(pr, rad(d), ps);
p[0] = ps[0] + pc[0];
p[1] = ps[1] + pc[1];
}
/* 誤差のチェック */
el = __abs(nl - l);
if (wflg == 1) d2 += 90;
p[2] = d2;
if (ep * 1e3 < el) { return -1; }
return 1;
}
void cimove(double c, double b, double a, int pm, double *p)
{
double x, y, r, p1, p2, w, d, q[2];
x = c[0];
y = c[1];
r = c[2];
p1 = c[3] ? c[3] : 0;
p2 = c[4] ? c[4] : 360;
w = c[5] ? c[5] : 1;
d = c[6] ? c[6] : 0;
if (pm == 0) {
w = (r * w + b) / (r + b);
r = r + b;
} else if (pm == -1 || pm == 1) {
if (pm == -1) a = -a;
polar(r, a, 1, q);
x += q[0]; y += q[1];
} else if (pm == 2) {
d = d + deg(a);
}
p[0] = x;
p[1] = y;
p[2] = r;
p[3] = p1;
p[4] = p2;
p[5] = w;
p[6] = d;
}
void line_x_poi(double *ln1, double *ln2, double *p)
{ /* 線と線の交点 */
double x1, y1, x2, y2, dx1, dy1, m1;
double x3, y3, x4, y4, dx2, dy2, m2;
double aa, b1, b2;
x1 = ln1[0];
y1 = ln1[1];
x2 = ln1[2];
y2 = ln1[3];
dx1 = x2 - x1;
dy1 = y2 - y1;
m1 = (dx1 == 0) ? 1e20 : dy1 / dx1;
x3 = ln2[0];
y3 = ln2[1];
x4 = ln2[2];
y4 = ln2[3];
dx2 = x4 - x3;
dy2 = y4 - y3;
m2 = (dx2 == 0) ? 1e20 : dy2 / dx2;
aa = m1 - m2;
if (__abs(aa) < 1e-10) aa = 0;
if (aa != 0) {
b1 = y1 - m1 * x1;
b2 = y3 - m2 * x3;
p[0] = (b2 - b1) / aa;
p[1] = (m1 * b2 - m2 * b1) / aa;
} else {
int flg = 0;
double p1[2], p2[2], p3[2], p4[2];
p1[0] = x1;
p1[1] = y1;
p2[0] = x2;
p2[1] = y2;
p3[0] = x3;
p3[1] = y3;
p4[0] = x4;
p4[1] = y4;
if (ptdist(p2, p3) < 1e-12) { p[0] = p2[0]; p[1] = p2[1]; flg = 1; }
if (ptdist(p1, p4) < 1e-12) { p[0] = p1[0]; p[1] = p1[1]; flg = 1; }
if (ptdist(p2, p4) < 1e-12) { p[0] = p2[0]; p[1] = p2[1]; flg = 1; }
if (ptdist(p1, p3) < 1e-12) { p[0] = p1[0]; p[1] = p1[1]; flg = 1; }
if (flg == 0) exit(NO_OPERATION);
}
}
double apa(int o, double p) {
return (90.0 * (2 * __int( (o + 1) / 2.0) - 1) - p) * pow(-1.0, o - 1.0);
}
double __fx__(double x, double k) {
return sqrt(1.0 - pow(k * sin(x), 2.0));
}
void gauss(int n) {
int i, j, e1;
double t, x0, pkm1, pk, t1, pkp1;
double den, d1, dpn, d2pn, d3pn, d4pn;
double u, v, h, pp, dp, x, f, w;
int m = __int( (n + n % 2) / 2);
e1 = n * (n + 1);
for (i = 1; i <= m; i++) {
t = (4.0 * i - 1.0) * PI / (4.0 * n + 2.0);
x0 = (1.0 - (1.0 - 1.0 / n) / (8.0 * n * n)) * cos(t);
pkm1= 1.0;
pk = x0;
for (j = 2; j <= n; j++) {
t1 = x0 * pk;
pkp1 = t1 - pkm1 - (t1 - pkm1) / j + t1;
pkm1 = pk;
pk = pkp1;
}
den = 1.0 - x0 * x0;
d1 = n * (pkm1 - x0 * pk);
dpn = d1 / den;
d2pn = (2.0 * x0 * dpn - e1 * pk) / den;
d3pn = (4.0 * x0 * d2pn + (2.0 - e1) * dpn) / den;
d4pn = (6.0 * x0 * d3pn + (6.0 - e1) * d2pn) / den;
u = pk / dpn;
v = d2pn / dpn;
h = -1.0 * u * (1.0 + 0.5 * u * (v + u * (v * v - d3pn / (3.0 * dpn))));
pp = pk + h * (dpn + 0.5 * h * (d2pn + h / 3.0 * (d3pn + 0.25 * h * d4pn)));
dp = dpn + h * (d2pn + 0.5 * h * (d3pn + h * d4pn / 3.0));
h -= pp / dp;
x = x0 + h;
f = d1 - h * e1 * (pk + 0.5 * h * (dpn + h / 3.0 * (d2pn + 0.25 * h * (d3pn + 0.2 * h * d4pn))));
w = 2.0 * (1.0 - x * x) / (f * f);
__ee__[i] = x;
__hh__[i] = w;
}
}
double elliptic_f(double p, double k, int n) {
double l = 0.0;
int i;
int m = __int( (n + n % 2) / 2);
if (k == 0 || n < 1) {
return p;
}
for (i = 1; i <= __int(n / 2.0); i++) {
l += __hh__[i] / __fx__(__ee__[i] * p, k);
}
if (m == __int(n / 2.0) + 1) { /* 奇数個の場合 */
l += __hh__[m] / __fx__(__ee__[m] * p, k) / 2.0;
}
return l * p;
}
double elliptic_e(double p, double k, int n) {
double l = 0.0;
int i;
int m = __int( (n + n % 2) / 2);
if (k == 0 || n < 1) {
return p;
}
for (i = 1; i <= __int(n / 2.0); i++) {
l += __hh__[i] * __fx__(__ee__[i] * p, k);
}
if (m == __int(n / 2.0) + 1) { /* 奇数個の場合 */
l += __hh__[m] * __fx__(__ee__[m] * p, k) / 2.0;
}
return l * p;
}
double arc_l(double r, double p1, double p2, double w) {
double p, lp;
int n, n1, n2, n3;
double k, ek;
p = p2 - p1;
while (p <= 0 ) p += 360;
while (p > 360) p -= 360;
if (w == 1) {
lp = p * PI / 180.0;
} else {
if (w > 1) {
w = 1.0 / w;
r /= w;
p1 = p1 - 90.0;
p2 = p2 - 90.0;
}
k = sqrt(1.0 - w * w);
n = __int(12 * (1.5 / sqrt(w))); /* ガウス積分点の数 */
if (n > 300) n = 300;
if (__abs(__ee__[0]) > 1 || __abs(__ee__[0]) < 0.5) gauss(n);
/* gauss(n); */
ek = elliptic_e(PI / 2.0, k, n);
lp = __int(__abs(p2 - p1) / 360) * 4.0 * ek;
while (p1 < 0 ) p1 += 360;
while (p1 >= 360) p1 -= 360;
while (p2 < 0 ) p2 += 360;
while (p2 >= 360) p2 -= 360;
if (angle_eql(p1, p2, 1e-12) == 0) {
return 4 * ek * r;
} else {
n1 = 1 + __int(p1 / 90.0);
if (n1 > 4) n1 = 4;
lp += ek * ( (n1 + 1) % 2) + elliptic_e(apa(n1, p1) * PI / 180, k, n) * pow(-1.0, n1 - 1.0);
n2 = 1 + __int(p2 / 90.0);
if (n2 > 4) n2 = 4;
if (n1 == n2) n1 += 4;
n3 = 5 + 2 * __int( (n2 - 1) / 2.0);
lp += ek * ( (n3 - n1) % 4) + elliptic_e(apa(n2, p2) * PI / 180, k, n) * pow(-1.0, (double) n2);
if (lp < 0.0) lp += 4.0 * ek;
if (lp > 4.0 * ek) lp -= 4.0 * ek;
}
}
return lp * r;
}
void arrow(FILE *f, double *pt1, double d, double size, double xang, double xan2)
{
double l, a, c, pt2[2], pt3[2], pt4[2];
l = size;
a = xang * PI / 180 / 2.0;
polarto(pt1, l, d - a, pt2);
polarto(pt1, l, d + a, pt4);
if (xan2 == 180) {
fprintf(f, "sl %.15g %.15g %.15g %.15g %.15g %.15g\n", pt1[0], pt1[1], pt2[0], pt2[1], pt4[0], pt4[1]);
} else {
c = cos(a) - sin(a) * tan( (90.0 - xan2 / 2.0) * PI / 180);
polarto(pt1, l * c, d, pt3);
fprintf(f, "sl %.15g %.15g %.15g %.15g %.15g %.15g %.15g %.15g\n", pt1[0], pt1[1], pt2[0], pt2[1], pt3[0], pt3[1], pt4[0], pt4[1]);
}
}
