// DY10
// 空気抵抗を受ける放射体の運動のグラフ
// 原作：平田邦男『新/BASICによる物理』共立出版
// 画面は原作(640×400)の2倍の1280×800。座標の数値も2倍にしてある。
// 原作の INPUT の代わりに、下の2つの値を書き換えて実行する。

// ---- 入力値（原作の INPUT に対応） ----
float v0 = 100;     // 初　速　度 [m/s]  // 空気抵抗を大きくすると初速を上げないといけない
float k = 0.005;   // 抵抗の係数 [1/m]  // 大きくすると空気抵抗が増えるY.O
// ------------------------------------

// 色（PC-98の色番号に対応）
color COL4 = color(0, 255, 0);    // 4: 緑
color COL5 = color(0, 255, 255);  // 5: 水色
color COL6 = color(255, 255, 0);  // 6: 黄

float pi = 3.1416, g = 9.8;
float x0 = 0, y0 = 0, dt = 0.01;
float x00 = 100, y00 = 700;  // 原点の座標
int d = 40, s = 10;          // d=目もり間隔, s=目もり線の長さ
int kk;
float o;  // 原作の STEP(2*S,O) の O（英字のオー）。値を入れていないので 0

void setup() {
  size(1280, 800);
  background(0);

  // 文字は1文字 16×32 として LOCATE の位置に置く
  textFont(createFont("SansSerif", 32));
  textAlign(LEFT, TOP);
  fill(COL5);
  // 原作で INPUT の問いと入力値が出ていた位置に、同じ形で表示する
  text("初　速　度 [m/s]：" + v0, 30*16, 4*32);
  text("抵抗の係数 [1/m]：" + k, 30*16, 5*32);

  stroke(COL4);
  line(x00, y00, x00 + 1100, y00);  // X軸
  for (int i = 1; i <= 27; i++) {   // 目もり付けのプログラム
    if (d*i % 200 == 0) line(x00 + d*i, y00, x00 + d*i, y00 - 2*s);
    line(x00 + d*i, y00, x00 + d*i, y00 - s);
  }
  line(x00, y00, x00, y00 - 600);   // Y軸
  for (int i = 1; i <= 14; i++) {
    if (d*i % 200 == 0) line(x00, y00 - d*i, x00 + 2*s, y00 - d*i + o);
    line(x00, y00 - d*i, x00 + s, y00 - d*i);
  }

  kk = 4;  // 座標軸の倍率
  text("一目もりの長さ = " + d/kk + " [m]", 30*16, 7*32);
  text("Y", 5*16, 4*32);
  text("0", 5*16, 22*32);
  text("X", 72*16, 22*32);

  stroke(COL6);
  strokeWeight(2);
  
  // 描きたい角度の前の//を外す
 // sub1(20);
 // sub1(30);
  sub1(45);
 // sub1(60);
 // sub1(70);
}

// サブルーチン
void sub1(float th) {
  float x = x0, y = y0;
  float vx = v0 * cos(th*pi/180), vy = v0 * sin(th*pi/180);
  float v = sqrt(vx*vx + vy*vy);
  float ax = -k*v*vx,             ay = -g - k*v*vy;
  vx = vx + ax*(dt/2);            vy = vy + ay*(dt/2);
  while (true) {
    x = x + vx*dt;                y = y + vy*dt;
    ax = -k*v*vx;                 ay = -g - k*v*vy;
    vx = vx + ax*dt;              vy = vy + ay*dt;
    v = sqrt(vx*vx + vy*vy);
    float x1 = kk*x + x00,        y1 = -kk*y + y00;
    if (x1 >= 1200 || y1 >= 700) break;
    point(x1, y1);
  }
}