「プログラムでシダを描画する」を SML# で シダ植物

夏休みの自由学習ということで ML をかじってみたわけですが、遅ればせながらプログラムを書いてみました。出来る限り元の Fortran プログラムの移植になるようにしました。

SML#には乱数や時計機能がないようなので、乱数による場合分けは、常に実行としました。

I/O がよく分からず、ものすごく遅いルーチンを書いてしまったので大きな図はかけませんでした。また再帰も遅いのでほどほどにで実行しています。

なお BMP ファイルの大きさは8の倍数である必要があります。

### 実行結果
main 12 240 240 "Shida.bmp"
![Shida_SML.png](https://qiita-image-store.s3.amazonaws.com/0/44522/e572e624-8a14-f754-a441-f04e58df35e4.png)

main 10 200 200 "Shida_10.bmp"
![Shida_ML_10.png](https://qiita-image-store.s3.amazonaws.com/0/44522/8e96f7d0-1929-69fd-34ae-b0f29482b2f8.png)

###プログラム
```　
(* BMP header I/O *)
fun bmp_file_header nx ny = 
  {bfType      = "BM", 
   bfSize      = 14 + 40 + nx * ny * 3, 
   bfReserved1 = 0, 
   bfReserved2 = 0, 
   bfOffBits   = 14 + 40}

fun bmp_info_header nx ny =
  {biSize          = 40, 
   biWidth         = nx , 
   biHeight        = ny, 
   biPlanes        = 1, 
   biBitCount      = 24, 
   biCompression   = 0, 
   biSizeImage     = nx * ny * 3, 
   biXPelsPerMeter = 3780, 
   biYPelsPerMeter = 3780, 
   biClrUsed       = 0, 
   biClrImportant  = 0}

fun str_int2 i = 
    str (chr (i mod 256)) 
  ^ str (chr (i div 256))  

fun str_int4 i = 
    str (chr (i                       mod 256))
  ^ str (chr (i div              256  mod 256))
  ^ str (chr (i div (      256 * 256) mod 256))
  ^ str (chr (i div (256 * 256 * 256) mod 256)) 

fun str_file_header fh =
    #bfType fh 
  ^ str_int4 (#bfSize      fh) 
  ^ str_int2 (#bfReserved1 fh) 
  ^ str_int2 (#bfReserved2 fh) 
  ^ str_int4 (#bfOffBits   fh)

fun str_info_header ih =
    str_int4 (#biSize          ih) 
  ^ str_int4 (#biWidth         ih)
  ^ str_int4 (#biHeight        ih)
  ^ str_int2 (#biPlanes        ih)
  ^ str_int2 (#biBitCount      ih)
  ^ str_int4 (#biCompression   ih)
  ^ str_int4 (#biSizeImage     ih)
  ^ str_int4 (#biXPelsPerMeter ih)
  ^ str_int4 (#biYPelsPerMeter ih)
  ^ str_int4 (#biClrUsed       ih)
  ^ str_int4 (#biClrImportant  ih)    

(* array I/O *)
fun str_array' a i n res =  (* DO LOOP :  DO i = _, n ... a(i) ... *)
  if i = n then res 
         else str_array' a (i + 1) n (res ^ (str (chr (Array.sub(a, i))))) 

fun str_array a =
  str_array' a 0 (Array.length a) ""

fun ipoint a nx ny ix iy (ir, ig, ib) =
let
  fun offset ix iy nx ny = 3 * ( ix + iy * nx ) 
  val ipos = offset ix iy nx ny
in 
  (
   Array.update(a, ipos    , ib); 
   Array.update(a, ipos + 1, ig); 
   Array.update(a, ipos + 2, ir)
  )
end;

(* plot point *)
fun point a nx ny x y = 
let 
  val ix = floor x  
  val iy = floor y  
in 
  ipoint a nx ny ix iy (0, 255, 0)
end;

(* fern *)
fun w1x x y =  0.836 * x + 0.044 * y;
fun w1y x y = ~0.044 * x + 0.836 * y + 0.169;
fun w2x x y = ~0.141 * x + 0.302 * y;
fun w2y x y =  0.302 * x + 0.141 * y + 0.127;
fun w3x x y =  0.141 * x - 0.302 * y;
fun w3y x y =  0.302 * x + 0.141 * y + 0.169;
fun w4x x y =  0.000 * x + 0.000 * y;
fun w4y x y =  0.000 * x + 0.175337 * y; 

fun fern' a nx ny i n x y =
  if i = n then point a nx ny ((Real.fromInt nx) * x * 0.98 + 0.5 * (Real.fromInt nx)) ((Real.fromInt ny) * y * 0.98)
           else (fern' a nx ny (i - 1) n (w1x x y) (w1y x y) ;
                 fern' a nx ny (i - 1) n (w2x x y) (w2y x y) ;
                 fern' a nx ny (i - 1) n (w3x x y) (w3y x y) ;
                 fern' a nx ny (i - 1) n (w4x x y) (w4y x y) )

(* main program *)
fun main n nx ny filename =
let
  val a = Array.array(nx * ny * 3, 0)
  val f = TextIO.openOut( filename )
in 
  (
   fern' a nx ny n 0 0.0 0.0;
   TextIO.output(f, str_file_header (bmp_file_header nx ny));
   TextIO.output(f, str_info_header (bmp_info_header nx ny));
   TextIO.output(f, str_array a);
   TextIO.closeOut(f)
  )
end;
```