FITS 转换为图像和图像透明度
将 FITS 文件转换为已知的图像格式，
引言
本文共三部分。
- 将FITS文件转换为更常见的图像格式
- 以动态透明度显示转换后的图像
- 播放延时视频，同时使用动态透明度。
必备组件
- FITS IO 库，可以在这里找到： http://fits.gsfc.nasa.gov/fits_libraries.html
- SDO FITS 文件，可以在这里找到： http://jsoc.stanford.edu/data/aia/synoptic/
- Visual Studio
背景
本文中的代码相当直接，但找到代码却有些困难。 在工作中，我们观察太阳，试图理解它（研究），并试图进行空间天气预测（操作）。
观察太阳的最佳工具之一是SDO卫星。（http://sdo.gsfc.nasa.gov/） 这颗卫星可以在不同的波段提供太阳日冕的观测（关于AIA仪器，还有其他仪器）。 AIA成像仪以几种波段呈现图像：094、131、171、193、211、304、335、1600、1700和4500埃。 每个波段显示不同的特征。 （取决于热量和密度）。
这个原型基于这样的想法：温度随着向外扩散而升高，穿过过渡区/日冕。
因此，通过根据基础温度（开尔文）对图像进行排序，我们应该可以得到某种可缩放的3D图像。
事实证明并非如此，但这个工具仍然很有用，因为图像都经过了相同的校准（对齐中心、对齐盘面大小等），因此我们可以更容易地通过“缩放”来跟踪不同波段的特征。
关于FITS
FITS代表Flexible Image Transport System（灵活图像传输系统）。 它是一种天文学中常用的标准文件格式。 它旨在处理表格和图像（表格也是一种图像）。 一个FITS文件可以包含多个表格和图像，因此文件的大小也会因来源而异。 目前C#对它的支持不是很好。 而且唯一的可用库不允许使用像idl或python库那样多的特殊功能。
使用代码
读取FITS
在C#中尝试实现“图像化”之前，您需要将FITS文件转换为.Net可以处理的格式：bmp、jpg、png、tiff等。
这相当困难，因为关于FITS文件包含什么信息并不多。 例如，每个仪器都可以提供这样的文件，但具有不同类型的值。 正如我通过SDO AIA FITS文件发现的那样，像素值不是RGB。
FITS文件包含HDU部分。 Header Data Units（头部数据单元）。 每个单元包含两个部分：一个包含原始数据的**数据**部分和一个包含元数据的**头部**部分。 一些**头部**关键字是强制的，其他是可选的（或可以自由添加）。
这是转换方法。 它接收FITS文件的文件路径。
（此转换特定于SDO AIA格式。）
public static System.Drawing.Image Convert(string pathfitsfile){
int totalmax = 0;
int totalmin = 0;
int wavelength = 0;
System.Drawing.Image result = null;
//http://fits.gsfc.nasa.gov/fits_libraries.html#CSharpFITS
//Read the fits file
nom.tam.fits.Fits fits = new nom.tam.fits.Fits(pathfitsfile);
nom.tam.fits.BasicHDU basichdu;
do{
basichdu = fits.readHDU();
if(basichdu != null){
basichdu.Info();
} //end if
} //end do
while(basichdu != null);
//loop through the Header Data Units
for(int i = 0; i < fits.NumberOfHDUs; i++){
basichdu = fits.getHDU(i);
//this retreives the current bandwidth
string card = basichdu.Header.GetCard(20);
wavelength = System.Convert.ToInt32(card.Replace("WAVELNTH=", "").Trim());
if(basichdu.GetType().FullName == "nom.tam.fits.ImageHDU"){
try{
nom.tam.fits.ImageHDU imghdu = (nom.tam.fits.ImageHDU)basichdu;
//get out the pixel array and the dimensions
Array [] array = (Array[])imghdu.Data.DataArray;
int x = imghdu.Axes[0];
int y = imghdu.Axes[1];
result = new System.Drawing.Bitmap(x, y);
System.Drawing.Graphics g = System.Drawing.Graphics.FromImage(result);
int idx_x = 0;
int idx_y = 0;
//get the real min/max values.
//We need it to get the values between 0 - 255.
for(int j = array.Length-1; j >= 0; j--){
int[] row = (int[])array[j];
for(int k = 0; k < row.Length; k++){
totalmax = (row[k] > totalmax ? row[k] : totalmax);
totalmin = (row[k] < totalmin ? row[k] : totalmin);
} //end for
}
//offset the values so everything is positive
int offset = (int)(totalmin < 0.0 ? (-1)*totalmin : 0.0);
//bottom is top, top is bottom.
for(int j = array.Length-1; j >= 0; j--){
//pixels are NOT RGB values so here we convert it to an RGB value.
idx_x = 0;
idx_y++;
int[] row = (int[])array[j];
for(int k = 0; k < row.Length; k++){
double val = row[k];
//This was trial and error.
switch(wavelength){
case 335:
case 94:
case 131:
break;
case 304: val /= 3;
break;
case 171:
case 193:
case 211:
case 1600: val /= 32;
break;
case 1700: val /= 128;
break;
case 4500: val /= 1024;
break;
default:
break;
} //end switch
//make sure all pixels are in range 0 - 255
val = (val > 255.0 ? val = 255.0 : val);
val = (val < 0.0 ? val = 0.0 : val);
System.Drawing.Color c = System.Drawing.Color.FromArgb(2013265920 + (int)val);
//call this function to color the image depending on wavelength
c = WeighRGBValue(System.Drawing.Color.FromArgb(c.A, c.B, c.B, c.B), wavelength);
((System.Drawing.Bitmap)result).SetPixel(idx_x, idx_y, c);
idx_x++;
}
}
} //end try
catch(Exception ex){
string error = ex.Message;
} //end catch
} //end if
} //end for
return result;
}
基本上，这段代码会剥离文件的**数据**部分（像素值），并将这些值转换为RGB。
如果您不使用 `WeighRGBValue` 方法，图像将是灰度的（这没问题）。 但为了后续缩放，您需要能够识别哪个特征属于哪个图像（或相应的）波段。 因此，我们对其进行着色。 我试图模拟网站上找到的原始SDO AIA颜色。
//http://www.rapidtables.com/web/color/RGB_Color.htm
private static System.Drawing.Color WeighRGBValue(System.Drawing.Color color, int wavelength){
//Tried to get a similar color to the original site.
System.Drawing.Color c;
int alpha, red, green, blue;
red = green = blue = 255;
alpha = color.A;
switch(wavelength){
case 94: //green : 00FF00
red = 0;
green = color.G;
blue = 0;
break;
case 131: //teal : 008080
red = 0;
green = color.G;
blue = color.B;
break;
case 304: //red : FF0000
red = color.R;
green = color.G / 5;
blue = color.B / 5;
break;
case 335: //blue : 0000FF
red = (color.R / 5);
green = (color.G / 5);
blue = color.B;
break;
case 171: //gold : FFD700
red = color.R;
green = (int)((double)color.G / 255.0 * 215.0);
blue = 0;
break;
case 193: //copper : B87333
red = color.R; //184 mapped to 255
green = (int)((double)color.G / 255.0 * 115.0); //115 mapped to 255
blue = (int)((double)color.B / 255.0 * 51.0); //50 mapped to 255
break;
case 211: //purple : 800080
red = color.R;
green = 0;
blue = color.B;
break;
case 1600: //ocher : BBBB00
red = (int)((double)color.R / 255.0 * 187.0); //BB mapped to 255
green = (int)((double)color.G / 255.0 * 187.0); //BB mapped to 255
blue = 0;
break;
case 1700: //pink : FFC0CB
red = color.R;
green = (int)((double)color.G / 255.0 * 192.0);
blue = (int)((double)color.B / 255.0 * 203.0);
break;
case 4500: //silver : C0C0C0
red = color.R;
green = color.G;
blue = color.B;
break;
default: red = 0;
green = 0;
blue = 0;
break;
} //end switch
c = System.Drawing.Color.FromArgb(alpha, red, green, blue);
return c;
}
简单来说，我们正在将灰色转移到另一种颜色，使其显示正确。 这涉及到一些试错。 我试图尽可能接近原始颜色，但也试图不丢失任何特征。
图像透明度
请注意，您也可以从SDO网站下载jpg图像，以跳过第一部分。
一旦您有了图像，就可以将它们一个一个地堆叠起来。 它们大小相同，太阳中心相同，盘面大小也相同。 我在项目源文件中添加了一行图像。
在WPF中，您为每个波段创建一个Image。 我使用温度作为基础值来排序我的图像。
- 304
- 131
- 171
- 193
- 211
- 335
- 1600
- 1700
- 4500
温度以开尔文（K）表示。
您可以在互联网上找到大部分信息，但这里也有一些
（不完全，但足以满足我们的需求）
| 波段 | K（从） | K（到） | Order |
|---|---|---|---|
| 0094 | 6300000 | 7 | |
| 0131 | 400000 | 16000000 | 2 |
| 0171 | 630000 | 3 | |
| 0193 | 1200000 | 20000000 | 4 |
| 0211 | 20000000 | 5 | |
| 0304 | 50000 | 1 | |
| 0335 | 25000000 | 6 | |
| 1600 | |||
| 1700 | |||
| 4500 |
排序基于K（从）。 值大部分是范围。（例如0131、0193）。
从磁盘加载图像可以这样做
string fn = rootfolder + @"\" + basefilename.Replace(bw, "0094");
images["0094"].BeginInit();
images["0094"].UriSource = new Uri(fn, UriKind.Absolute);
images["0094"].EndInit();
img_0094.Source = images["0094"];
img_0094.Stretch = Stretch.Uniform;
这段代码加载0094图像。 对所有波段重复此操作。 我在XAML中手动添加和样式化了Image控件，但您也可以轻松地在代码中添加，并使用List或Array来保存控件。
请注意，每个波段的文件名是相似的（因此可以重用文件名并插入带有新波段的波段部分）。 您还可以使用文件名来更改“时间”参数，稍后在“延时视频”部分中使用。
为了缩放，我添加了一个滑块，并在其 `ValueChanged` 事件中添加了此代码。
//0 - 1000
if(slider_zoom.Value >= 0 && slider_zoom.Value < 100){
img_0304.Opacity = (100 - slider_zoom.Value)/100;
img_0131.Opacity = (0 + slider_zoom.Value) / 100;
}
if(slider_zoom.Value >= 100 && slider_zoom.Value < 200){
img_0131.Opacity = (100 - (slider_zoom.Value - 100))/100;
img_0171.Opacity = (0 + (slider_zoom.Value - 100))/100;
}
if(slider_zoom.Value >= 200 && slider_zoom.Value < 300){
img_0171.Opacity = (100 - (slider_zoom.Value - 200))/100;
img_0193.Opacity = (0 + (slider_zoom.Value - 200))/100;
}
if(slider_zoom.Value >= 300 && slider_zoom.Value < 400){
img_0193.Opacity = (100 - (slider_zoom.Value - 300))/100;
img_0211.Opacity = (0 + (slider_zoom.Value - 300))/100;
}
if(slider_zoom.Value >= 400 && slider_zoom.Value < 500){
img_0211.Opacity = (100 - (slider_zoom.Value - 400))/100;
img_0335.Opacity = (0 + (slider_zoom.Value - 400))/100;
}
if(slider_zoom.Value >= 500 && slider_zoom.Value < 600){
img_0335.Opacity = (100 - (slider_zoom.Value - 500))/100;
img_0094.Opacity = (0 + (slider_zoom.Value - 500))/100;
}
if(slider_zoom.Value >= 600 && slider_zoom.Value < 700){
img_0094.Opacity = (100 - (slider_zoom.Value - 600))/100;
img_1600.Opacity = (0 + (slider_zoom.Value - 600))/100;
}
if(slider_zoom.Value >= 700 && slider_zoom.Value < 800){
img_1600.Opacity = (100 - (slider_zoom.Value - 700))/100;
img_1700.Opacity = (0 + (slider_zoom.Value - 700))/100;
}
if(slider_zoom.Value >= 800 && slider_zoom.Value < 900){
img_1700.Opacity = (100 - (slider_zoom.Value - 800))/100;
img_4500.Opacity = (0 + (slider_zoom.Value - 800))/100;
}
if(slider_zoom.Value >= 900 && slider_zoom.Value < 1000){
img_4500.Opacity = (100 - (slider_zoom.Value - 900))/100;
}
我必须承认，这比我想象的要容易，而且进行得很顺利。
延时视频
这部分仍在建设中。 当然，任何想法都欢迎。
这是我的想法
基本上，每个波段需要两条轨道（双缓冲），在后台加载下一个图像，同时显示前一个。 这意味着您要处理20个Image控件，10个显示，10个在后台。（可能的内存问题）
我认为问题在于在循环遍历图像时动态更改 `Opacity` 属性。 由于这绝对是多线程的，所以可能会成为一个问题。
不用说，如果我能让它工作，我会在本节中发布。
关注点
这个原型项目的有趣之处在于它非常直观。 处理图像总是这样的。 您稍微改变一下代码，就能立即看到结果。
这是一个“空闲时间”项目，所以我也不确定能花多少时间，但我会尽力完成它。
关于fits转换
这被证明是一个很大的挑战。 C#可用的FITS库不是最好的，而且我了解到并非所有图像都由RGB像素组成。 将其转换为所需的格式被证明是一个巨大的挑战。 最终结果“尚可”，但仍然不像在SDO网页上看到的那样。 我非常确定我的试错过程可以改进，特别是跳过（或设置）阈值以上的像素值（在转换之前！）为黑色。 我看到一些python代码里是这么做的。 效果是太阳周围的阴霾应该会有些消失。 我们用来除以FITS像素值的任意值也可以改进。 我猜FITS头部中可能有一些缩放因子可以使用，但到目前为止我还没有找到。
关于图像透明度
起初，我以为你需要对图像本身进行操作才能实现这个功能，而不是改变控件的透明度。 由于这并不需要太多工作，我只是试了一下，效果非常好而且流畅。 不用说，这也可以用于其他类型的图像，在这些图像中，做类似的事情可能很有用。
延时视频
创建延时视频本身不应该太难。 您创建一个计时器，然后一个接一个地加载/显示图像。 在运行延时视频时尝试控制图像可能会是另一回事，但让10个波段同步运行同时玩转透明度将是一个巨大的挑战。
参考文献
- “由NASA/SDO以及AIA、EVE和HMI科学团队提供。”
历史
- 2016-04-08：版本1.0