# coding: utf-8
# The aim of this code is to demonstrate how one might apply an ocean colour algorithm to level 2 reflectances in order to estimate chlorophyll-a concentration.

# Import a number of required modules and functions to aid in the processing
import shutil, os, sys
import argparse
import numpy as np
import matplotlib.pyplot as plt
import netCDF4 as nc



# No we define the input data from level 2 (the example used here shows eddy features in the chl field but you could change this to any other level 2 data file set for a region of your own interest).
# We copy the file to the local /tmp space to avoid overwriting the orignal data.
#example area of the Yellow Sea
indir='./DATA/S3A_OL_2_WFR____20181002T020947_20181002T021247_20181003T112306_0179_036_217_2340_MAR_O_NT_002.SEN3'
COPY_DIR = indir.replace('.SEN3', '_copy.SEN3')
COPY_DIR= './tmp/'+os.path.basename(COPY_DIR)
if os.path.isdir(COPY_DIR):
    print(COPY_DIR+' already exists')
else:
    shutil.copytree(indir, COPY_DIR)

# Here we load the required netcdfs
Oa06_reflectance=nc.Dataset(COPY_DIR+'/Oa06_reflectance.nc', 'r')
Oa06_reflectance_values=nc.Dataset(COPY_DIR+'/Oa06_reflectance.nc', 'r').variables['Oa06_reflectance']
Oa07_reflectance_values=nc.Dataset(COPY_DIR+'/Oa07_reflectance.nc', 'r').variables['Oa07_reflectance']
Oa08_reflectance_values=nc.Dataset(COPY_DIR+'/Oa08_reflectance.nc', 'r').variables['Oa08_reflectance']
Oa10_reflectance_values=nc.Dataset(COPY_DIR+'/Oa10_reflectance.nc', 'r').variables['Oa10_reflectance']

# #Shen 2010 equations 2,4 and 5
Hchl=(0.74*Oa10_reflectance_values[:] + 0.26*Oa07_reflectance_values[:] - Oa08_reflectance_values[:])

HDelta=Oa07_reflectance_values[:]-0.5*(Oa06_reflectance_values[:]+Oa10_reflectance_values[:])

SCI=Hchl-HDelta

# # Now we convert SCI to Chl
CHL_SCI_Spring=179378.0*pow(SCI,2) +92.934*SCI + 0.2736
# # We will write the file out the new band to a file
# print('about to write')

dims=[x for x in Oa06_reflectance.dimensions]
ds = nc.Dataset('./tmp/test_output.nc','w', format="NETCDF4")
for dimension in Oa06_reflectance.dimensions:
    ds.createDimension(dimension)
ds.createVariable('Chl_SCI_spring', 'f8',dimensions=dims, fill_value=CHL_SCI_Spring.get_fill_value())
ds.variables['Chl_SCI_spring'][:]=CHL_SCI_Spring[:]
ds.close()
# print('finished writing')