HiTRANで透過法のガススペクトルを計算しCSVで出力させる

# 背景

フーリエ変換型赤外分光光度計（FTIR）を使用して有機ガスの赤外スペクトルを測定した際、ガスの典型的な振動回転スペクトルが得られる。
塩化水素の振動回転スペクトルは学生実験でも採用しやすい典型的な振動回転スペクトルが得られる。(図１)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3673555/c623be9f-c061-449d-8301-1ad0e9f8f3ea.png)


:::note info
上記データは日本分光のFT/IR-6100、光路長が10cmのガラス製のガスセル、窓板がNaClで測定したデータ　...だったと思います。当時の実験ノートを見ても窓板は流石に記録してなかったです。
分解能は0.5cm-1、積算回数は32回で測定しています。データ見れば分かりますがDTGS検出器です。
:::


同様に他の有機ガスで回転振動スペクトルを得ようとする。通常はガスセルを使うのだが、その窓板,光路長は種類が大量にあり、 **予め選択しなければならない。** 

選択をするためにはモル吸光係数εと濃度が必要になる。
モル吸光係数εは量子化学計算を行えばある程度推定することができる。
しかし、量子化学計算は計算リソースが必要となり、ソフトによっては有料になる。

:::note info
ガスの振動回転スペクトルの計算ぐらいであればそこまでリソースがいらないので、gaussianのライセンス料が払えるかが重要。　
GAMESSなどを使えば無料でも出来るが情報が少なめな点には注意。
:::

量子化学計算を使わない方法として市販の定量用ガスライブラリを使う方法もある。
濃度と光路長の情報がライブラリに含まれるためモル吸光係数εを算出するやり方である。
これは有料ではあるが上記よりかは安価に行える。

しかしライブラリを使った方法ではアサイメントが難しくなる。そこでアサイメントとスペクトルの計算の両方ができるデータベース「HiTRAN」を使った方法を記述する。

:::note warn
HiTRANそのものは多くの論文によるデータベースではあるものの、アサイメントを含む論文をpublishする際はGaussian等の量子化学計算をリファレンスにしたものにすることを推奨します。
:::

# 計算方法　

:::note alert
以降の作業は各種公式の[リファレンスマニュアル](https://hitran.org/static/hapi/hapi_manual.pdf)を読むべきです。
コードはあくまで一例です。
Windowsでの作業を想定しています。Macでも原理的には動作するはずですが未検証です。
:::

### Python環境構築 

[Python公式ページ](https://www.python.org/downloads/windows/ "Python Releases for Windows")からPythonのインストーラーをダウンロードし、インストール作業を行う。
インストールの際にパスを通すことを忘れずに。

:::note info
Pythonは人気の言語なのでQiitaで検索すればより詳しいインストール方法が学べるはずです。HiTRANではPythonの実行環境が必須です。

インストールが嫌ならばGoogle Colaboratoryの使用をオススメいたします。
:::

### パッケージのインストール

[HiTRAN公式ページ](https://hitran.org/hapi/ "HiTRAN HAPI download")から `hapi.py`をダウンロードし、pythonのコードを置いている,置くつもりのディレクトリに設置する。

...もしくは
Google Colaboratoryなどローカル以外の環境でPythonを動かしている場合は以下の`pip`コマンドを使用してインストールする。 Google Colabの場合は`!pip install hitran-api`と打つ

```
pip install hitran-api
```

### importする
HAPIは以下のようにimportすれば使える。

```python
from hapi import *
```

HiTRANのデータベースへの接続は以下の文で行う。

```python
db_begin('.')
```
データベースへ接続後、`fetch`関数に任意の引数を入れて呼び出すことで必要なデータをダウンロードすることができる。

```python
fetch('HCl', 14, 1, 400, 4000)
```

<details><summary>引数</summary>

第一引数　`データ名`　
上記の例では取得したデータはHCl.dataに格納される。

第二引数　`moleculeID`
HiTRANにおける分子のID　[Line-by-Line Search](https://hitran.org/lbl/)のページでIDを確認できる。

第三引数　`isotopologueID`
HiTRANにおける各分子の同位体のID　目的の分子にチェックを入れたあと`Select Isotopologues`をクリックしたあとに表示されるページから確認可能　[HClの例](https://hitran.org/lbl/2?15=on) 

</details>

上記の例では原子量が35のHClのデータをダウンロードしている。（存在比75.76%）
実際のスペクトルと比べるときは、原子量35のHClと原子量37のHClのコンボリューションスペクトルを計算して比較するのが望ましい。

`fetch`関数を使用してデータを取得した場合はローカルにデータが残るので、同じ分子であれば２回目以降は`fetch`関数を呼び出さなくても良い。

### データを処理する

得られた`hogehoge.data`は多量のパラメーターを持ったデータ群である。このうち必要なものを抜き出してからデータ処理する必要がある。
今回は最終が波数対吸光度のグラフを得たいので、波数と吸光度を計算するためのパラメーターが必要になる。
ガスのスペクトルを計算するには、

+ 圧力
+ 温度
+ 濃度

があれば吸収断面積（モル吸光係数）が求められる。

求めるには`absorptionCoefficient_Lorentz()`関数を呼び出す

```python
absorptionCoefficient_Lorentz(SourceTables='HCl', 
                                         HITRAN_units=False, 
                                          WavenumberStep=0.1, 
                                          Environment={'p': 1.0, 'T': 298}, 
                                          Diluent={'air': 1.0-c,'HCl':c})
```

<details><summary>引数</summary>

引数`SourceTables`　利用する`データ名`　
上記の例では取得したデータHCl.dataを利用する

引数`HITRAN_units` 論理型　`TRUE` or `FALSE`
TRUE： `cm^2/molec`
FALSE: `cm-1`
吸収スペクトルや透過率などの計算はFALSEに設定する

引数`WavenumberStep` 数値型
デフォルトは0.01cm-1
あとで`SciPy`ライブラリで分解能のリサンプリングを行う方針であれば、0.001以下に設定する
（そもそもそのやり方の方がより厳密な計算が行える。）

引数`Environment` 
第一引数　圧力[atm]　デフォルト＝1.0
第二引数　温度[K]　デフォルト＝296.0

引数`Diluent` 混合ガスの定義
A+B+C...=1になるように定義する

</details>

### 装置関数を使用（やりたければ）

装置関数を使用することで実際のスペクトルに近づけることができる。例えば装置のアポダイゼーションが三角関数である場合は以下のような関数を実行してデータを変換する。

```python
nu,coef,i1,i2,Res = convolveSpectrum(nu, coef, SlitFunction=SLIT_TRIANGULAR,Resolution=1.0)
```


### 吸光度を計算する

吸光度は`Lambert-Beer則`に基づいて計算する。下記は10.0cmガスセルを使用して濃度が100ppmとした時の計算

```python
absorbance = coef * 10.0 * 1.0e-4
```
あとは`numpy`ライブラリ、`matplotlib`ライブラリ、`pandas`ライブラリうを使っていい感じにデータを表示して出力させる。

# 具体例 塩化水素のスペクトル

```python:HiTRAN_HCl_csv.py
# 必要なパッケージをインポート
import numpy as np
from hapi import *
import matplotlib.pyplot as plt
import pandas as pd

# HITRANデータベースへの接続を開始
db_begin('.')

# 塩化水素（HCl）のデータをダウンロード
fetch('HCl', 15, 1, 400, 4000)

# パラメータを設定
P = 1.0 # 圧力 (atm)
T = 298.0 # 温度 (K)
L = 10.0 # 光路長 (cm)
c = 1.0e-4 # 濃度 (ppm)

# 吸収係数の計算
nu, coef = absorptionCoefficient_Lorentz(SourceTables='HCl', 
                                         HITRAN_units=False, 
                                          WavenumberStep=0.1, 
                                          Environment={'p': P, 'T': T}, 
                                          Diluent={'air': 1.0-c,'HCl':c})

# 装置関数の使用
nu,coef,i1,i2,Res = convolveSpectrum(nu, coef, SlitFunction=SLIT_TRIANGULAR,Resolution=1.0)

#吸光度の計算
absorbance = coef * L * c

# グラフを作成
plt.figure(figsize=(10, 6))
plt.plot(nu, absorbance)

# 軸のラベルを設定
plt.xlabel('wavelength (cm-1)')
plt.ylabel('Absorbance (Abs.)')

# グラフのタイトルを設定
plt.title('Absorbance spectrum of HCl')

# グラフを表示
plt.show()

# データフレームを作成
df = pd.DataFrame({
    'wavelength': nu,
    'Absorbance': absorbance
})

# CSVファイルに保存
df.to_csv('spectrum_HCl.csv', index=False)

```

得られた計算スペクトルは以下の通りです。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3673555/23b53045-8508-43d9-a0eb-02c3b274bd7c.png)

csvファイルも同じフォルダに出力されてていい感じです。

一方で最初の実験スペクトルとのAbs.解離に気が付いた方もいるかもしれません。これは実験でのHClガスの捕集方法として、ガスセルを真空引きして塩酸試薬瓶に近づけて捕集する方法を使用したので、おそらく100ppm以上あると思います。今回は適当に100ppmで計算をしましたが、他の分析方法例えばGC-MSやTOF-MSを使用して濃度をはっきりさせてから計算スペクトルを計算させるとかなり近い計算スペクトルが得られると考えられます