PythonでFITSファイルを効率的にフィルタリングする方法（不等号条件にも対応）

はじめに

FITSファイルは天文学で広く使われるデータフォーマットであり、観測データを管理・解析するために非常に便利です。Pythonでは、astropyライブラリを使って簡単にFITSファイルを読み書きし、データのフィルタリングも可能です。今回は、不等号条件（>=, <, <=, >）にも対応したフィルタリング方法を紹介します。

初心者でも分かりやすいように、擬似データを使ってFITSファイルを生成し、そのデータに対して条件付きフィルタを適用する具体例を見ていきます。

本記事の実行例は下記の google Colab にあげています。

---

FITSファイルを用いた例

1. 擬似データの作成とFITSファイルの生成

まず、擬似的なデータセットを作成し、それをFITSファイルに保存します。astropy.io.fitsを使うことで、簡単にFITSファイルを作成できます。

from astropy.io import fits
import numpy as np

# 擬似データの作成
data = np.zeros(100, dtype=[('COLUMN_A', 'f4'), ('COLUMN_B', 'f4')])
data['COLUMN_A'] = np.random.uniform(0, 100, size=100)
data['COLUMN_B'] = np.random.uniform(0, 500, size=100)

# FITSファイルの生成
hdu = fits.BinTableHDU(data)
hdu.writeto('sample_data.fits', overwrite=True)
print("FITSファイル 'sample_data.fits' を生成しました。")

解説:

• np.zerosで100行のデータを作成し、COLUMN_Aには0〜100の乱数を、COLUMN_Bには0〜500の乱数をそれぞれ生成しました。
• fits.BinTableHDUを使って、バイナリテーブル形式のFITSファイルを生成し、ファイルに保存しています。

これで、FITSファイル sample_data.fits が作成されました。

---

2. フィルタ条件のパース

次に、カラムに対するフィルタ条件を指定します。例えば、「COLUMN_Aが50以上で、COLUMN_Bが300未満」という条件を指定します。これらの条件を扱いやすい形式にパースするための関数を作成します。

import re

def parse_filter_conditions(conditions):
    filters = []
    # 正規表現でカラム、演算子、値を分離    
    # '!=' を含めた正規表現パターン
    condition_pattern = re.compile(r"(.*?)(==|<=|>=|<|>|!=)(.*)")
    for condition in conditions.split(","):
        match = condition_pattern.match(condition.strip())
        if match:
            col, op, value = match.groups()
            filters.append((col.strip(), op, float(value.strip())))
        else:
            raise ValueError(f"Invalid filter condition: {condition}")
    return filters

解説:

• この関数は、フィルタ条件を「カラム名」「演算子」「値」に分割します。例えば、COLUMN_A>=50 という条件は、COLUMN_A（カラム名）、>=（演算子）、50（値）に分割されます。
• 正規表現を使うことで、複数のフィルタ条件をカンマ区切りで指定可能です。

---

3. フィルタの適用

続いて、フィルタ条件を実際にデータに適用します。条件に基づいてデータを絞り込みます。

# フィルタ適用関数に '!=' の処理を追加
def apply_filters(data, filters):
    mask = np.ones(len(data), dtype=bool)
    for col, op, value in filters:
        if op == "==":
            mask &= (data[col] == value)
        elif op == "!=":
            mask &= (data[col] != value)
        elif op == "<":
            mask &= (data[col] < value)
        elif op == "<=":
            mask &= (data[col] <= value)
        elif op == ">":
            mask &= (data[col] > value)
        elif op == ">=":
            mask &= (data[col] >= value)
    return data[mask]    

解説:

• maskは、データの各行が条件を満たしているかどうかを表すブール配列です。
• 各フィルタ条件に応じてmaskを更新し、最終的に条件を満たす行のみを返します。

---

4. FITSファイルの読み込みとフィルタリングの実行

最後に、作成したFITSファイルを読み込み、指定した条件に基づいてフィルタリングを行います。

# FITSファイルの読み込みとフィルタ適用
with fits.open('sample_data.fits') as hdul:
    data = hdul[1].data  # テーブルデータを取得
    print("FITSファイルのデータを読み込みました。")
    
    # フィルタ条件を設定
    conditions = "COLUMN_A>=50,COLUMN_B<300"
    filters = parse_filter_conditions(conditions)
    
    # フィルタを適用
    filtered_data = apply_filters(data, filters)
    
    # 結果の表示
    print("フィルタ適用後のデータ:")
    for row in filtered_data:
        print(row)

解説:

• fits.openでFITSファイルを読み込み、データを取得します。
• フィルタ条件として「COLUMN_Aが50以上かつCOLUMN_Bが300未満」という条件を指定し、フィルタリングを実行します。
• フィルタ適用後のデータが表示されます。

---

辞書に対して適用する場合

FITSファイルの場合に、apply_filters関数が正しく動作する理由は、FITSファイルのデータ構造がNumPyの構造化配列（Structured Array）であるためです。

FITSファイルのデータ構造

FITSファイルのデータは、astropy.io.fitsモジュールを使って読み込むと、NumPyの構造化配列として扱われます。構造化配列は、カラム名（列名）を使ってデータにアクセスできるだけでなく、各カラムがNumPy配列として保存されています。これにより、以下のようにカラムごとに配列を扱う形になります。

data['COLUMN_A']  # これで 'COLUMN_A' カラムの NumPy 配列が取得できる

構造化配列は、FITSファイル内のテーブル形式のデータを扱うために非常に適しています。apply_filters関数がFITSデータに対して動作するのは、この理由です。data[col]の部分がNumPy配列を返すため、NumPyのブール演算（==, <, <=, >, >=）が直接利用できます。

辞書の場合

一方で、辞書形式の場合は、各カラムがリストや配列として保存されている可能性がありますが、辞書全体にはNumPyの構造化配列のようなデータ構造はありません。そのため、辞書で同様の操作を行おうとすると、各カラムのサイズが一致していなかったり、データ構造に違いが生じて、エラーが発生します。

辞書はFITSファイルの構造化配列と異なるため、辞書形式のデータに対しては別途フィルタリングの処理が必要です。下記に辞書に適用する場合のコードの例を挙げておきます。

サンプルコード

import numpy as np
import re

def parse_filter_conditions(conditions):
    filters = []
    condition_pattern = re.compile(r"(.*?)(==|<=|>=|<|>|!=)(.*)")
    for condition in conditions.split(","):
        match = condition_pattern.match(condition.strip())
        if match:
            col, op, value = match.groups()
            filters.append((col.strip(), op, float(value.strip())))
        else:
            raise ValueError(f"Invalid filter condition: {condition}")
    return filters

def apply_filters_dic(data, filters):
    mask = np.ones(len(list(data.values())[0]), dtype=bool)  # 各カラムのサイズに基づくマスクを作成
    for col, op, value in filters:
        if op == "==":
            mask &= (data[col] == value)
        elif op == "!=":
            mask &= (data[col] != value)
        elif op == "<":
            mask &= (data[col] < value)
        elif op == "<=":
            mask &= (data[col] <= value)
        elif op == ">":
            mask &= (data[col] > value)
        elif op == ">=":
            mask &= (data[col] >= value)
    return {col: values[mask] for col, values in data.items()}

# ダミーデータの作成
data = {
    'COLUMN_A': np.array([10, 20, 30, 40, 50]),
    'COLUMN_B': np.array([100, 200, 300, 400, 500]),
}

# フィルタ条件 ('!=' を含む例)
conditions = "COLUMN_A>=20,COLUMN_B!=300"

# 条件をパース
filters = parse_filter_conditions(conditions)

# フィルタを適用
filtered_data = apply_filters_dic(data, filters)

# 結果の表示
print("Filtered Data:")
for key, value in filtered_data.items():
    print(f"{key}: {value}")

まとめ

今回の解説では、Pythonのastropyライブラリを使って、FITSファイルに対するフィルタリングを行う方法を紹介しました。不等号条件にも対応しているため、複雑なフィルタリングが可能です。

• フィルタ条件の指定：不等号を使った柔軟なフィルタ条件の指定が可能です。
• フィルタの適用：複数のフィルタ条件に基づいてデータを絞り込むことができます。

この方法を使えば、効率的にFITSファイル内のデータを解析・抽出できます。