# -*- coding: utf-8 -*-
"""
Created on Thu Jan 14 18:00:32 2021
@author: Zdenek Navratil

Sample program for simulation IR spectra using Hitran database
 """

import numpy as np
import matplotlib.pyplot as mp
from hapi import *

"""****************************************************************************
HAPI download - neni nutné 
****************************************************************************"""
def DownloadAll():
  for i in range(1, 50, 1):
      fetch ( '{:02d}_01'.format(i) , i, 1, 500, 4000) # h_{:02d}_01


'''---------------------------------------------------------------------------
 Výpis dostupných molekul a popis dat 
---------------------------------------------------------------------------'''
getHelp(ISO_ID)  

#molec_id -- číslo molekuly (je v databázi čar)
#local_iso_id (číslo izotopu dané molekuly podle zastoupení (abundance)
#ISO_ID je řazeno dle "global" isotopologue number

# vyber molekuly
myID=5

#moleculeName(mol_id)
print(moleculeName(myID))

#abundance(mol_id, local_iso_id)
print(abundance(myID,1))  # mol_id, iso_id
print(abundance(myID,2))  # 

#isotopologueName((mol_id, local_iso_id)
print(isotopologueName(myID, 2)) 


print(ISO_ID[26])

'''---------------------------------------------------------------------------
 Vytvoření lokální databáze
---------------------------------------------------------------------------'''

db_begin('dbx')  

#fetch(TableName, M, I, numin, numax)
#fetch_by_ids(TableName, iso_id_list, numin, numax)

fetch('db_1', myID , 1, 500, 4000)
fetch('db_2', myID , 2, 500, 4000)
#fetch_by_ids('CO-p2', [26], 500, 4000)

describe('db_1')

X, Y = getColumns ('db_1',['nu', 'sw'])

mp.figure()
mp.title("Database integral cross section data")
mp.xlabel("wavenumber  (cm$^{-1}$)")
mp.ylabel('S')
mp.xlim(2000, 2300)
mp.plot(X, Y)
mp.show()

tableList()

'''---------------------------------------------------------------------------
Spočtení absorpčního koeficientu 
---------------------------------------------------------------------------'''
#absorptionCoefficient_Voigt
#absorptionCoefficient_Lorentz
#absorptionCoefficient_Doppler
#absorptionCoefficient_HT  #  using HT (Hartmann Tran) profile. https://www.sciencedirect.com/science/article/pii/S0022407313002422?via%3Dihub


wn1, xs1 = absorptionCoefficient_Voigt([(myID,1)], ['db_1'], Environment = {'p':1.0,'T':296.0}, Diluent={'air':1, 'self':1e-5}, OmegaRange=(500,4000), OmegaStep=0.01, HITRAN_units=False  )
wn2, xs2 = absorptionCoefficient_Voigt([(myID,1),(myID,2)], ['db_1', 'db_2'], Environment = {'p':1.0,'T':296.0}, Diluent={'air':1, 'self':1e-5}, OmegaRange=(500,4000), OmegaStep=0.01, HITRAN_units=False  )

fig, ax = mp.subplots()
mp.xlabel("wavenumber  (cm$^{-1}$)")
mp.title("Absorption coefficient")
mp.plot(wn2, xs2, label='1 + 2')
mp.plot(wn2, 10*(xs2-xs1), label='10 x difference')
mp.xlim(2000, 2300)
ax.invert_xaxis()
mp.legend()
mp.show()

#Partiční sumy
partitionSum(myID,1,296)
partitionSum(myID,2,296)

'''---------------------------------------------------------------------------
Spočtení spektra
---------------------------------------------------------------------------'''

concentration = 1e-4 # relative to 1 bar

wn3, a3 =  transmittanceSpectrum(wn1, concentration*xs1, Environment={'l':20.0, 'T':296.0})
# l is length of the cuvette in cm

sps=[(wn3, a3, 0)]

for res in [0.5,2.0,4.0]:
  wni, ai, i1, i2, slit =  convolveSpectrum(wn3, a3, SlitFunction=SLIT_GAUSSIAN, Resolution=res)
  sps = sps + [(wni.copy(), ai.copy(), res)]
  


fig, ax = mp.subplots()
mp.xlabel("wavenumber  (cm$^{-1}$)")
mp.ylabel("transmittance")
mp.title("Transmittance spectrum at different instrumental resolution")
for s in sps:
        mp.plot(s[0], s[1], label='{:.1f} {:s}'.format(s[2], "cm$^{-1}$"))

mp.xlim(2175, 2185)



ax.invert_xaxis()
mp.legend() 
mp.show()