ピクセルを操作してマンデルブロ集合の画像を作る

#はじめに
画像処理でピクセルを操作する練習として、マンデルブロ集合(詳しくは[Wikipedia](https://en.wikipedia.org/wiki/Mandelbrot_set)で)の画像を作ってみた。Javaは全体像、C#とC++は適当に拡大した部分の画像ができる。
基本的には、`image.setPixel(i, j, color)`のようなメソッドを使うよりも、内部の配列に直接書き込む方が速いらしい。

#Java
C#やC++がチャンネルごとの配列になっているのに対して、JavaはARGBを一つの32ビット整数にまとめたものが要素になっている。

```java
import java.awt.Color;
import java.awt.image.BufferedImage;
import javax.imageio.ImageIO;
import java.io.File;
import java.io.IOException;

public class Mandelbrot{
    public static final int WIDTH = 1000;
    public static final int HEIGHT = 1000;
        
    public static void main(String[] args){
		int[] pixels = new int[WIDTH*HEIGHT];
		for(int i = 0; i < WIDTH; i++){
			for(int j = 0; j < HEIGHT; j++){
				pixels[WIDTH*j+i] = calcMandelbrot(i, j);
			}
		}
		BufferedImage image = new BufferedImage(WIDTH, HEIGHT, BufferedImage.TYPE_INT_ARGB);
		image.setRGB(0, 0, WIDTH, HEIGHT, pixels, 0, WIDTH);
		File file = new File("mandelbrot.png");
		try{
			ImageIO.write(image, "png", file);
		}catch(IOException e){
			System.out.println(e);
		}
	}
	
	// マンデルブロ集合のパラメータ
    public static final double xMin = -2.1;
    public static final double xMax = 0.5;
    public static final double yMin = -1.3;
    public static final double yMax = 1.3;
    public static final int nMax = 300;
     
    // 点(i, j)に対する色をAHSVで計算
	public static int calcMandelbrot(int i, int j){
		final double c = xMin+i*(xMax-xMin)/WIDTH;
		final double d = yMin+j*(yMax-yMin)/HEIGHT;
        double x1 = 0.0, y1 = 0.0, x2, y2;
        int n = 0;
        for(n = 0; n < nMax; n++){
            x2 = x1*x1-y1*y1+c;
            y2 = 2*x1*y1+d;
            if(x2*x2+y2*y2 > 4.0) break;
            x1 = x2;
            y1 = y2;
        }
		final float t = (float)n/nMax;
		if(t >= 1.0f){
			return 0xff000000; //black
		}else{
			return 0xff000000 | Color.HSBtoRGB(t, 0.6f, 1.0f);
		}
	}
}
```

結果：
![mandelbrot1.png](https://qiita-image-store.s3.amazonaws.com/0/161876/1e3770bf-81a9-6e26-b881-34acec6867fb.png)


参考：[Javaで画像の画素を操作する | なべろぐ](https://www.google.co.jp/url?sa=t&rct=j&q=&esrc=s&source=web&cd=1&cad=rja&uact=8&ved=0ahUKEwi5s-2q-4XSAhWKgrwKHYoDAugQFggaMAA&url=http%3A%2F%2Fblog.nabe-ch.net%2F%3Fp%3D532&usg=AFQjCNHl-RIE4VzT1sbyfrpzRNMPOLTUeA&sig2=ckRLdb2H-aPoCA82xN8r2g)

#C# 
C#にはHSVとRGBとの変換がないらしいので、[Wikipedia](https://ja.wikipedia.org/wiki/HSV%E8%89%B2%E7%A9%BA%E9%96%93)を参考にして作った。

```csharp
using System;
using System.Drawing;
using System.Drawing.Imaging;

namespace Mandelbrot{
    static class Program{
        public const int WIDTH = 1000;
        public const int HEIGHT = 1000;

        public static void Main(string[] args){
            Bitmap image = new Bitmap(WIDTH, HEIGHT, PixelFormat.Format32bppArgb);
            BitmapData data = image.LockBits(
                new Rectangle(0, 0, WIDTH, HEIGHT),
                ImageLockMode.WriteOnly,
                PixelFormat.Format32bppArgb);
            byte[] buf = new byte[4*WIDTH*HEIGHT];
            for(int i = 0; i < WIDTH; i++){
                for(int j = 0; j < HEIGHT; j++){
                    byte a = 0, r = 0, g = 0, b = 0;
                    CalcMandelbrot(i, j, ref a, ref r, ref g, ref b);
                    buf[4*(WIDTH*j+i)+3] = a;
                    buf[4*(WIDTH*j+i)+2] = r;
                    buf[4*(WIDTH*j+i)+1] = g;
                    buf[4*(WIDTH*j+i)] = b;
                }
            }
            System.Runtime.InteropServices.Marshal.Copy(buf, 0, data.Scan0, buf.Length);
            image.UnlockBits(data);
            image.Save("mandelbrot.png");
        }
        
        // マンデルブロ集合のパラメータ
        public const double xMin = -0.7406219098542647;
        public const double xMax = -0.7406219098519411;
        public const double yMin = 0.15805475052205210;
        public const double yMax = 0.15805475052421664;
        public const int nMax = 18000;
        
        // 点(i, j)に対する色をAHSVで計算
        public static void CalcMandelbrot(int i, int j, ref byte a, ref byte r, ref byte g, ref byte b){
            double c = xMin+i*(xMax-xMin)/WIDTH;
            double d = yMin+j*(yMax-yMin)/HEIGHT;
            double x1 = 0.0, y1 = 0.0, x2, y2;
            int n = 0;
            for(n = 0; n < nMax; n++){
                x2 = x1*x1-y1*y1+c;
                y2 = 2*x1*y1+d;
                if(x2*x2+y2*y2 > 4.0) break;
                x1 = x2;
                y1 = y2;
            }
            double t = (double)n/nMax;
            if(t >= 1.0){
                a = 0xff;
                r = g = b = 0; //black
            }else{
                ARGBfromAHSV(ref a, ref r, ref g, ref b, 1.0, 0.55, 0.3+0.3*Math.Sin(12*Math.PI*t), 0.7+0.3*Math.Cos(16*Math.PI*t));
            }
        }
        
        // ARGBからAHSVへ変換
        public static void ARGBfromAHSV(ref byte aOut, ref byte r, ref byte g, ref byte b, double aIn, double h, double s, double v){
            aOut = (byte)(255*aIn);
            r = (byte)(255*v);
            g = (byte)(255*v);
            b = (byte)(255*v);
            if (s <= 0.0) return;
            h *= 6.0;
            int i = (int)h;
            double f = h-i;
            switch(i){
                case 0:
                    g = (byte)(g*(1-s*(1-f)));
                    b = (byte)(b*(1-s));
                    break;
                case 1:
                    r = (byte)(r*(1-s*f));
                    b = (byte)(b*(1-s));
                    break;
                case 2:
                    r = (byte)(r*(1-s));
                    b = (byte)(b*(1-s*(1-f)));
                    break;
                case 3:
                    r = (byte)(r*(1-s));
                    g = (byte)(g*(1-s*f));
                    break;
                case 4:
                    r = (byte)(r*(1-s*(1-f)));
                    g = (byte)(g*(1-s));
                    break;
                case 5:
                    g = (byte)(g*(1-s));
                    b = (byte)(b*(1-s*f));
                    break;
            }
        }
    }
}
```

結果：
![mandelbrot2.png](https://qiita-image-store.s3.amazonaws.com/0/161876/9db70542-7282-0afe-1cdf-da6c59ac8d86.png)


参考：[[C#] ビットマップにピクセル単位で高速にアクセスするには (GetPixel/SetPixel vs BitmapData 速度比較)](http://www.84kure.com/blog/2014/07/13/c-%E3%83%93%E3%83%83%E3%83%88%E3%83%9E%E3%83%83%E3%83%97%E3%81%AB%E3%83%94%E3%82%AF%E3%82%BB%E3%83%AB%E5%8D%98%E4%BD%8D%E3%81%A7%E9%AB%98%E9%80%9F%E3%81%AB%E3%82%A2%E3%82%AF%E3%82%BB%E3%82%B9/)

#C++
OpenCVを使った。今回は全部不透明にしてるので特に意味ないけど、OpenCVではalpha値を含めた場合のHSVとRGBとの変換ができないらしい([c++ - OpenCV cv::cvtColor drops alpha channel, How to keep alpha data?](http://stackoverflow.com/questions/20229577/opencv-cvcvtcolor-drops-alpha-channel-how-to-keep-alpha-data))。

```cpp
#include <cmath>
#include <opencv2/opencv.hpp>
#include <opencv2/imgproc.hpp>

void calcMandelbrot(int i, int j, unsigned char& h, unsigned char& s, unsigned char& v);

constexpr int width = 1000;
constexpr int height = 1000;

int main(){
	cv::Mat image(height, width, CV_8UC3);
	for(int i = 0; i < width; i++){
		for(int j = 0; j < height; j++){
			unsigned char h, s, v;
			calcMandelbrot(i, j, h, s, v);
			image.data[3*(width*j+i)] = h;
			image.data[3*(width*j+i)+1] = s;
			image.data[3*(width*j+i)+2] = v;
		}
	}
	cv::cvtColor(image, image, CV_HSV2RGB);
	cv::imwrite("mandelbrot.png", image);
	return 0;
}

// マンデルブロ集合のパラメータ
constexpr double xMin = -1.4467108885212991;
constexpr double xMax = -1.4467096751794892;
constexpr double yMin = -1.3212921606783162e-5;
constexpr double yMax = -1.2078843332240347e-5;
constexpr int nMax = 1000;

// 点(i, j)に対する色をHSVで計算
// hの範囲は0から180まで
void calcMandelbrot(int i, int j, unsigned char& h, unsigned char& s, unsigned char& v){
	const double c = xMin+i*(xMax-xMin)/width;
	const double d = yMin+j*(yMax-yMin)/height;
	double x1 = 0.0, y1 = 0.0, x2, y2;
	int n = 0;
	for(n = 0; n < nMax; n++){
		x2 = x1*x1-y1*y1+c;
		y2 = 2*x1*y1+d;
		if(x2*x2+y2*y2 > 4.0) break;
		x1 = x2;
		y1 = y2;
	}
	const double t = static_cast<double>(n)/nMax;
	if(t >= 1.0){
		h = s = v = 0; //black
	}else{
		h = static_cast<unsigned char>(180*0.3);
		s = static_cast<unsigned char>(255*(0.3+0.3*sin(12*M_PI*t)));
		v = static_cast<unsigned char>(255*(0.7+0.3*cos(16*M_PI*t)));
	}
}
```

結果：
![mandelbrot3.png](https://qiita-image-store.s3.amazonaws.com/0/161876/cb9538fc-3d57-2ae9-ff03-7f1b81f6d04a.png)


参考：
・[Learning OpenCV 2](http://www.wakayama-u.ac.jp/~chen/opencv/opencv2.html)
・[緑色のピクセルを数える](http://www.cellstat.net/color/)