/*---------------------------------------------------------------------------*/
/* gccT2CNATM-gnuplot.c (2心円NATMトンネルモデル) */
/*---------------------------------------------------------------------------*/
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>
#define pi M_PI
 
void WGPL(char fnameW[50],char fnameF[50],int ksp,double B,double H,double W,double Ti,double ts,double tp,double Ae,double Ac,double Ls);
/*---------------------------------------------------------------------------*/
int main(int argc,char *argv[])
{
    char fnameW[50];
    char fnameF[50];
    double B,H,W,Ti,ts,tp,phi;
    double rr,hh,dw,RL,theta,A1,A2,A3,Ae,Ac,Ls;
    int ksp;

    strcpy(fnameW,argv[1]); /*コマンド記述ファイル名*/
    strcpy(fnameF,argv[2]); /*出力画像ファイル名*/

    ksp=atoi(argv[3]);  /*-1:モデル図，0:NATMトンネル*/
    B  =atof(argv[4]);  /*内空幅*/
    H  =atof(argv[5]);  /*内空高*/
    W  =atof(argv[6]);  /*インバート内空幅*/
    Ti =atof(argv[7]);  /*インバート厚*/
    ts =atof(argv[8]);  /*吹付厚*/
    tp =atof(argv[9]);  /*設計線から支払線までの距離*/

    /* 掘削面積 (Pay line)*/
        rr=0.5*B;
        hh=H-rr;
        dw=rr-0.5*W;
        RL=0.5*(hh*hh+dw*dw)/dw;
 
        rr=0.5*B+ts+tp;
        hh=H-rr+Ti+tp;
        theta=asin(hh/(RL+ts+tp));
        A1=(RL+ts+tp)*(RL+ts+tp)*(theta-sin(theta)*cos(theta));
        A2=(B+2.0*ts+2.0*tp-(RL+ts+tp)*(1.0-cos(theta)))*hh;
        A3=0.5*pi*rr*rr;
        Ae=A1+A2+A3;
        Ae=Ae*1e-6;
    /* インバートコン面積 (Concrete)*/
        rr=0.5*B;
        hh=H-rr;
        theta=asin(hh/RL);
        Ac=0.5*(W+W-2.0*Ti*tan(theta))*Ti;
        Ac=Ac*1e-6;
    /* 吹付長 (Shotcrete)*/
        rr=0.5*B+ts;
        hh=H-rr+Ti;
        theta=asin(hh/(RL+ts));
        Ls=rr*pi+2.0*(RL+ts)*theta;
        Ls=Ls*1e-3;

    printf("Ae=%.1f Ac=%.1f Ls=%.1f\n",Ae,Ac,Ls);
    
    WGPL(fnameW,fnameF,ksp,B,H,W,Ti,ts,tp,Ae,Ac,Ls);
    return 0;
}
/*---------------------------------------------------------------------------*/
void WGPL(char fnameW[50],char fnameF[50],int ksp,double B,double H,double W,double Ti,double ts,double tp,double Ae,double Ac,double Ls)
{
    FILE *fp;
    double xx,yy,x1,y1,x2,y2,x3,y3,x4,y4,xc,yc,phi;
    double BB,HU,HL,ds,dss,dl;
    double rr,hh,dw,RL,theta;
    char str[50],sv[50];
    int i,n=20;
 
    BB=B+2.0*ts+2.0*tp;
    HU=H+ts+tp;
    HL=(Ti+tp);
    ds=1000.0;
    dss=200.0;
    dl=400.0;
 
    rr=0.5*B;
    hh=H-rr;
    dw=rr-0.5*W;
    RL=0.5*(hh*hh+dw*dw)/dw;
 
    fp=fopen(fnameW,"w");
    fprintf(fp,"reset\n");
    fprintf(fp,"set terminal postscript eps enhanced font \"Helvetica\" 16\n");
    fprintf(fp,"set output \"%s\"\n",fnameF);
    fprintf(fp,"#\n");
    fprintf(fp,"set multiplot\n");
    fprintf(fp,"#\n");
    fprintf(fp,"set origin 0,0\n");
    fprintf(fp,"set size ratio -1\n");
    fprintf(fp,"#\n");
    fprintf(fp,"set xrange [-7000:7000]\n");
    fprintf(fp,"set yrange [-4000:10000]\n");
    fprintf(fp,"set border 0\n");
    fprintf(fp,"set xtics 0,0,0\n");
    fprintf(fp,"set ytics 0,0,0\n");
    fprintf(fp,"set xtics nomirror\n");
    fprintf(fp,"set ytics nomirror\n");
    fprintf(fp,"#\n");
    fprintf(fp,"set style arrow 1 size graph 0.02,10 filled nohead linewidth 1 linetype 5\n");/*一点鎖線*/
    fprintf(fp,"set style arrow 2 size graph 0.02,10 filled nohead linewidth 1 linetype 1\n");/*実線*/
    fprintf(fp,"set style arrow 3 size graph 0.02,10 filled heads linewidth 1 linetype 1\n"); /*両矢印*/
    fprintf(fp,"set style arrow 4 size graph 0.02,10 filled linewidth 1 linetype 1\n");       /*片矢印*/
    /*一点鎖線*/
    x1= 0.5*BB+ds+dss ; y1=H-0.5*B;
    x2=-0.5*BB-0.5*ds ; y2=y1;
    x3=0.0            ; y3=-HL-0.5*ds;
    x4=0.0            ; y4= HU+0.5*ds;
    fprintf(fp,"set arrow as 1 from %.3f,%.3f to %.3f,%.3f\n",x1,y1,x2,y2);
    fprintf(fp,"set arrow as 1 from %.3f,%.3f to %.3f,%.3f\n",x3,y3,x4,y4);
    fprintf(fp,"#\n");
    /*実線*/
        x1=0.0           ;y1=H+ts;
        x2=0.5*BB+ds+dss ;y2=y1;
    fprintf(fp,"set arrow as 2 from %.3f,%.3f to %.3f,%.3f\n",x1,y1,x2,y2);
        x1=0.0           ;y1=H;
        x2=0.5*BB+2.0*ds+dss ;y2=y1;
    fprintf(fp,"set arrow as 2 from %.3f,%.3f to %.3f,%.3f\n",x1,y1,x2,y2);
        x1=0.0           ;y1=0.0;
        x2=0.5*BB+2.0*ds+dss ;y2=y1;
    fprintf(fp,"set arrow as 2 from %.3f,%.3f to %.3f,%.3f\n",x1,y1,x2,y2);
        x1=0.0           ;y1=-Ti;
        x2=0.5*BB+ds+dss ;y2=y1;
    fprintf(fp,"set arrow as 2 from %.3f,%.3f to %.3f,%.3f\n",x1,y1,x2,y2);
        x1=-0.5*B  ;y1=H-0.5*B;
        x2=x1      ;y2=-HL-2.0*ds-dss;
    fprintf(fp,"set arrow as 2 from %.3f,%.3f to %.3f,%.3f\n",x1,y1,x2,y2);
        x1=0.5*B  ;y1=H-0.5*B;
        x2=x1     ;y2=-HL-2.0*ds-dss;
    fprintf(fp,"set arrow as 2 from %.3f,%.3f to %.3f,%.3f\n",x1,y1,x2,y2);
        x1=-0.5*W  ;y1=0.0;
        x2=x1      ;y2=-HL-ds-dss;
    fprintf(fp,"set arrow as 2 from %.3f,%.3f to %.3f,%.3f\n",x1,y1,x2,y2);
        x1=0.5*W   ;y1=0.0;
        x2=x1      ;y2=-HL-ds-dss;
    fprintf(fp,"set arrow as 2 from %.3f,%.3f to %.3f,%.3f\n",x1,y1,x2,y2);
 
    /*両矢印*/
        sprintf(sv,"%.0f",0.5*B);
        strcpy(str,sv);
        strcpy(sv," rotate by 90");
        x1=0.5*BB+ds  ;y1=H-0.5*B;
        x2=x1         ;y2=H;
        x3=x1-dl      ;y3=0.5*(y1+y2);
    fprintf(fp,"set arrow as 3 from %.3f,%.3f to %.3f,%.3f\n",x1,y1,x2,y2);
    if(ksp<0){strcpy(str,"r");strcpy(sv,"");}
    fprintf(fp,"set label \"%s\" at %.3f,%.3f center%s\n",str,x3,y3,sv);
        sprintf(sv,"%.0f",hh);
        strcpy(str,sv);
        strcpy(sv," rotate by 90");
        x1=0.5*BB+ds  ;y1=0.0;
        x2=x1         ;y2=H-0.5*B;
        x3=x1-dl      ;y3=0.5*(y1+y2);
    fprintf(fp,"set arrow as 3 from %.3f,%.3f to %.3f,%.3f\n",x1,y1,x2,y2);
    if(ksp<0){strcpy(str,"h");strcpy(sv,"");}
    fprintf(fp,"set label \"%s\" at %.3f,%.3f center%s\n",str,x3,y3,sv);
        sprintf(sv,"%.0f",H);
        strcpy(str,sv);
        strcpy(sv," rotate by 90");
        x1=0.5*BB+2.0*ds ;y1=0.0;
        x2=x1            ;y2=H;
        x3=x1-dl         ;y3=0.5*(y1+y2);
    fprintf(fp,"set arrow as 3 from %.3f,%.3f to %.3f,%.3f\n",x1,y1,x2,y2);
    if(ksp<0){strcpy(str,"H");strcpy(sv,"");}
    fprintf(fp,"set label \"%s\" at %.3f,%.3f center%s\n",str,x3,y3,sv);
        sprintf(sv,"%.0f",W);
        strcpy(str,sv);
        x1=-0.5*W ; y1=-HL-ds;
        x2= 0.5*W ; y2=y1;
        x3=0.5*(x1+x2);y3=y2-dl;
    fprintf(fp,"set arrow as 3 from %.3f,%.3f to %.3f,%.3f\n",x1,y1,x2,y2);
    if(ksp<0)strcpy(str,"W");
    fprintf(fp,"set label \"%s\" at %.3f,%.3f center\n",str,x3,y3);
        sprintf(sv,"%.0f",B);
        strcpy(str,sv);
        x1=-0.5*B     ; y1=-HL-2.0*ds;
        x2= 0.5*B     ; y2=y1;
        x3=0.5*(x1+x2); y3=y2-dl;
    fprintf(fp,"set arrow as 3 from %.3f,%.3f to %.3f,%.3f\n",x1,y1,x2,y2);
    if(ksp<0)strcpy(str,"B");
    fprintf(fp,"set label \"%s\" at %.3f,%.3f center\n",str,x3,y3);
 
    /*片矢印*/
        sprintf(sv,"%.0f",ts);
        strcpy(str,sv);
        strcpy(sv," rotate by 90");
        x1=0.5*BB+ds  ;y1=H+ts+0.7*ds;
        x2=x1         ;y2=H+ts;
        x3=x1+dl      ;y3=0.5*(y2+H);
    fprintf(fp,"set arrow as 4 from %.3f,%.3f to %.3f,%.3f\n",x1,y1,x2,y2);
    if(ksp<0){strcpy(str,"t{/=14 s}");strcpy(sv,"");}
    fprintf(fp,"set label \"%s\" at %.3f,%.3f center%s\n",str,x3,y3,sv);
        sprintf(sv,"%.0f",Ti);
        strcpy(str,sv);
        strcpy(sv," rotate by 90");
        x1=0.5*BB+ds  ;y1=-(Ti+0.7*ds);
        x2=x1         ;y2=-Ti;
        x3=x1+dl      ;y3=0.5*(y2);
    fprintf(fp,"set arrow as 4 from %.3f,%.3f to %.3f,%.3f\n",x1,y1,x2,y2);
    if(ksp<0){strcpy(str,"T{/=14 I}");strcpy(sv,"");}
    fprintf(fp,"set label \"%s\" at %.3f,%.3f center%s\n",str,x3,y3,sv);
        sprintf(sv,"%.0f",dw);
        strcpy(str,sv);
        strcpy(sv," rotate by 90");
        x1=-0.5*B-0.7*ds ;y1=-HL-ds;
        x2=-0.5*B        ;y2=y1;
        x3=0.5*(x2-0.5*W);y3=y1-dl;
    fprintf(fp,"set arrow as 4 from %.3f,%.3f to %.3f,%.3f\n",x1,y1,x2,y2);
    if(ksp<0)strcpy(str,"dw");
    fprintf(fp,"set label \"%s\" at %.3f,%.3f center\n",str,x3,y3);
        sprintf(sv,"%.0f",dw);
        strcpy(str,sv);
        x1=0.5*B+0.7*ds  ;y1=-HL-ds;
        x2=0.5*B         ;y2=y1;
        x3=0.5*(x2+0.5*W);y3=y1-dl;
    fprintf(fp,"set arrow as 4 from %.3f,%.3f to %.3f,%.3f\n",x1,y1,x2,y2);
    if(ksp<0)strcpy(str,"dw");
    fprintf(fp,"set label \"%s\" at %.3f,%.3f center\n",str,x3,y3);
 
        strcpy(str,"r=");
        sprintf(sv,"%.0f",0.5*B);
        strcat(str,sv);
        x1=0.0;
        y1=H-0.5*B;
        x2=rr*cos(45.0/180.0*pi);
        y2=y1+rr*sin(45.0/180.0*pi);
        x3=0.5*x2-dl*sin(45.0/180.0*pi);
        y3=y1+0.5*rr*sin(45.0/180.0*pi)+dl*cos(45.0/180.0*pi);
    fprintf(fp,"set arrow as 4 from %.3f,%.3f to %.3f,%.3f\n",x1,y1,x2,y2);
    if(ksp<0)strcpy(str,"r");
    fprintf(fp,"set label \"%s\" at %.3f,%.3f center rotate by 45\n",str,x3,y3);
 
        strcpy(str,"R=");
        sprintf(sv,"%.0f",RL);
        strcat(str,sv);
        theta=0.5*asin(hh/RL)/pi*180.0;
        x1=(RL-rr)-(RL-0.45*W)*cos(theta/180.0*pi);
        y1=(H-0.5*B)-(RL-0.45*W)*sin(theta/180.0*pi);
        x2=(RL-rr)-RL*cos(theta/180.0*pi);
        y2=(H-0.5*B)-RL*sin(theta/180.0*pi);
        x3=(RL-rr)-(RL-0.1*W)*cos(theta/180.0*pi)-dl*sin(theta/180.0*pi);
        y3=(H-0.5*B)-(RL-0.1*W)*sin(theta/180.0*pi)+dl*cos(theta/180.0*pi);
    fprintf(fp,"set arrow as 4 from %.3f,%.3f to %.3f,%.3f\n",x1,y1,x2,y2);
    if(ksp<0)strcpy(str,"R");
    fprintf(fp,"set label \"%s\" at %.3f,%.3f left rotate by %.3f\n",str,x3,y3,theta);
 
    if(ksp<0){
            theta=45.0;
            x1=       -(0.5*B+ts-ds)*cos(theta/180.0*pi);
            y1=H-0.5*B+(0.5*B+ts-ds)*sin(theta/180.0*pi);
            x2=       -(0.5*B+ts)*cos(theta/180.0*pi);
            y2=H-0.5*B+(0.5*B+ts)*sin(theta/180.0*pi);
        fprintf(fp,"set arrow as 4 from %.3f,%.3f to %.3f,%.3f\n",x1,y1,x2,y2);
            strcpy(str,"t{/=14 p}");
            x1=       -(0.5*B+ts+tp+ds)*cos(theta/180.0*pi);
            y1=H-0.5*B+(0.5*B+ts+tp+ds)*sin(theta/180.0*pi);
            x2=       -(0.5*B+ts+tp)*cos(theta/180.0*pi);
            y2=H-0.5*B+(0.5*B+ts+tp)*sin(theta/180.0*pi);
            x3=0.5*(x1+x2)+dl*sin(theta/180.0*pi);
            y3=0.5*(y1+y2)+dl*cos(theta/180.0*pi);
        fprintf(fp,"set arrow as 4 from %.3f,%.3f to %.3f,%.3f\n",x1,y1,x2,y2);
        fprintf(fp,"set label \"%s\" at %.3f,%.3f center rotate by %.3f\n",str,x3,y3,-theta);
    }
    fprintf(fp,"#\n");
 
    if(0<=ksp){
        /*数量描画*/
        x3=0.0;
            sprintf(sv,"%.1f",Ae);
            strcpy(str,"Ae=");
            strcat(str,sv);
            strcat(str," m^{/=10 2}");
            y3=8000.0;
        fprintf(fp,"set label \"%s\" at %.3f,%.3f\n",str,x3,y3);
            sprintf(sv,"%.1f",Ac);
            strcpy(str,"Ac=");
            strcat(str,sv);
            strcat(str," m^{/=10 2}");
            y3=7500.0;
        fprintf(fp,"set label \"%s\" at %.3f,%.3f\n",str,x3,y3);
            sprintf(sv,"%.1f",Ls);
            strcpy(str,"Ls=");
            strcat(str,sv);
            strcat(str," m");
            y3=7000.0;
        fprintf(fp,"set label \"%s\" at %.3f,%.3f\n",str,x3,y3);
    }
    fprintf(fp,"#\n");
 
    /*構造線描画*/
    fprintf(fp,"plot \"-\" with lines linewidth 5 linetype 1 notitle,\\\n");      /*吹付外形線*/
    if(ksp<0)fprintf(fp,"\"-\" with lines linewidth 3 linetype 2 notitle,\\\n"); /*支払い線*/
    fprintf(fp,"\"-\" with lines linewidth 5 linetype 1 notitle,\\\n");           /*吹付内形線*/
    fprintf(fp,"\"-\" with lines linewidth 5 linetype 1 notitle,\\\n");           /*吹付内形線*/
    fprintf(fp,"\"-\" with lines linewidth 5 linetype 1 notitle\n");              /*インバート線*/
    /*吹付外形(一挙に書く)*/
        xc=RL-rr;
        yc=H-0.5*B;
        yy=hh+Ti;
        theta=asin(yy/(RL+ts))/pi*180.0;
        phi=0.0;
        for(i=0;i<=n;i++){
            phi=theta/(double)n*(double)i;
            xx=xc-(RL+ts)*cos(phi/180.0*pi);
            yy=yc-(RL+ts)*sin(phi/180.0*pi);
            fprintf(fp,"%.3f %.3f\n",xx,yy);
        }
        xc=rr-RL;
        yc=H-0.5*B;
        phi=0.0;
        for(i=n;i>=0;i--){
            phi=theta/(double)n*(double)i;
            xx=xc+(RL+ts)*cos(phi/180.0*pi);
            yy=yc-(RL+ts)*sin(phi/180.0*pi);
            fprintf(fp,"%.3f %.3f\n",xx,yy);
        }
        fprintf(fp,"e\n");
 
    if(ksp<0){
        /*支払い線(一挙に書く)*/
        xc=RL-rr;
        yc=H-0.5*B;
        yy=hh+Ti+tp;
        theta=asin(yy/(RL+ts+tp))/pi*180.0;
        phi=0.0;
        for(i=0;i<=n;i++){
            phi=theta/(double)n*(double)i;
            xx=xc-(RL+ts+tp)*cos(phi/180.0*pi);
            yy=yc-(RL+ts+tp)*sin(phi/180.0*pi);
            fprintf(fp,"%.3f %.3f\n",xx,yy);
        }
        xc=rr-RL;
        yc=H-0.5*B;
        phi=0.0;
        for(i=n;i>=0;i--){
            phi=theta/(double)n*(double)i;
            xx=xc+(RL+ts+tp)*cos(phi/180.0*pi);
            yy=yc-(RL+ts+tp)*sin(phi/180.0*pi);
            fprintf(fp,"%.3f %.3f\n",xx,yy);
        }
        fprintf(fp,"e\n");
    }
 
    /*吹付内形*/
        xc=RL-rr;
        yc=H-0.5*B;
        yy=hh+Ti;
        theta=asin(yy/RL)/pi*180.0;
        phi=0.0;
        for(i=0;i<=n;i++){
            phi=theta/(double)n*(double)i;
            xx=xc-RL*cos(phi/180.0*pi);
            yy=yc-RL*sin(phi/180.0*pi);
            fprintf(fp,"%.3f %.3f\n",xx,yy);
        }
        fprintf(fp,"e\n");
        xc=rr-RL;
        yc=H-0.5*B;
        phi=0.0;
        for(i=0;i<=n;i++){
            phi=theta/(double)n*(double)i;
            xx=xc+RL*cos(phi/180.0*pi);
            yy=yc-RL*sin(phi/180.0*pi);
            fprintf(fp,"%.3f %.3f\n",xx,yy);
        }
        fprintf(fp,"e\n");
    /*インバート線*/
        xx=-0.5*W; yy=0.0;
        fprintf(fp,"%.3f %.3f\n",xx,yy);
        xx= 0.5*W; yy=0.0;
        fprintf(fp,"%.3f %.3f\n",xx,yy);
        fprintf(fp,"e\n");
 
    /*円形部描画*/
    fprintf(fp,"set origin 0,0\n");
    fprintf(fp,"set parametric\n");
    fprintf(fp,"set trange [0:pi]\n");
    fprintf(fp,"plot %.3f*cos(t),%.3f+%.3f*sin(t) with lines linewidth 5 linetype 1 notitle,\\\n",rr,H-0.5*B,rr);    /*内空*/
    if(ksp<0)fprintf(fp,"     %.3f*cos(t),%.3f+%.3f*sin(t) with lines linewidth 3 linetype 2 notitle,\\\n",rr+ts+tp,H-0.5*B,rr+ts+tp); /*支払い線*/
    fprintf(fp,"     %.3f*cos(t),%.3f+%.3f*sin(t) with lines linewidth 5 linetype 1 notitle\n",rr+ts,H-0.5*B,rr+ts); /*掘削径*/
    fprintf(fp,"#\n");
    fprintf(fp,"#\n");
    fprintf(fp,"unset multiplot\n");
    fprintf(fp,"#\n");
    fprintf(fp,"quit\n");
    fprintf(fp,"#\n");
    fclose(fp);
 
}
/*---------------------------------------------------------------------------*/