Commit 9a10d52e authored by CEVAER's avatar CEVAER
Browse files

Added variable attributes

parent ee6d0940
......@@ -2,10 +2,12 @@ import xarray as xr
import numpy as np
from tc_center.polar import to_polar, to_north_ref
from tc_center.plotting import contour_indices
from tc_center.functions import CenterStatus
import owi
from shapely.geometry import LineString, Point
import os
import logging
logger = logging.getLogger(__name__)
......@@ -46,22 +48,37 @@ def save_netcdf(sat_file, path, center_status, percent_outside, percent_inside_i
ds_polar_north["eye_center"] = xr.DataArray(data=Point(lon_center, lat_center).wkt,
attrs={"long_name": "TC eye center"})
ds_polar_north["low_wind_area_barycenter"] = xr.DataArray(data=Point(lon_eye_lwind, lat_eye_lwind).wkt,
attrs={"long_name": "Low wind area barycenter"})
attrs={"long_name": "Low wind area barycenter (WKT)"})
ds_polar_north["high_wind_area_barycenter"] = xr.DataArray(data=Point(lon_eye_hwind, lat_eye_hwind).wkt,
attrs={"long_name": "High wind area barycenter"})
attrs={"long_name": "High wind area barycenter (WKT)"})
ds_polar_north["interpolated_track_point"] = xr.DataArray(data=Point(lon_eye_fg, lat_eye_fg).wkt,
attrs={"long_name": "Closest interpolated track point"})
ds_polar_north["track_point_status"] = xr.DataArray(data=",".join([str(ct_st.value) for ct_st in center_status]))
ds_polar_north["percent_outside"] = xr.DataArray(data=percent_outside)
ds_polar_north["percent_inside_island"] = xr.DataArray(data=percent_inside_island)
ds_polar_north["percent_non_usable"] = xr.DataArray(data=percent_non_usable)
#ds_polar_north.attrs["TC center status"] = status
#ds_polar_north.attrs["Area within 100km around track point outside acquisition bounding box in percent"] = percent_outside
#ds_polar_north.attrs["Area of 100km around track point inside acquisition bounding and on land in percent"] = percent_inside_island
#ds_polar_north.attrs["% of area of 100km around track point inside " \
attrs={"long_name": "Closest interpolated track point"
" (WKT)"})
ds_polar_north["track_point_flag"] = xr.DataArray(data=",".join([str(ct_st.value) for ct_st in center_status]),
attrs={"long_name": "Track point flags",
"flag_values": "inside, outside, land",
"flag_meanings": "inside: track point is inside "
"satellite acquisition. outside : track"
" point is outside acquisition. land: "
"track point is on land."})
ds_polar_north["percent_outside"] = xr.DataArray(data=percent_outside,
attrs={"long_name": "Area within 100km around track point outside"
" acquisition bounding box, in percent"})
ds_polar_north["percent_inside_island"] = xr.DataArray(data=percent_inside_island,
attrs={"long_name": "Area of 100km around track point inside"
" acquisition bounding box and on land, "
"in percent"})
ds_polar_north["percent_non_usable"] = xr.DataArray(data=percent_non_usable, attrs={
"long_name": "Sum of percent_outside and percent_inside_island, which is the area of 100km around track point"
"which is unusable for analysis, in percent."
})
# ds_polar_north.attrs["TC center status"] = status
# ds_polar_north.attrs["Area within 100km around track point outside acquisition bounding box in percent"] = percent_outside
# ds_polar_north.attrs["Area of 100km around track point inside acquisition bounding and on land in percent"] = percent_inside_island
# ds_polar_north.attrs["% of area of 100km around track point inside " \
# "acquisition bounding and on land or outside of bounding box"] = percent_inside_island
lon_lwind, lat_lwind = research_contour(msk_lwind, lon_algo, lat_algo)
......@@ -69,14 +86,18 @@ def save_netcdf(sat_file, path, center_status, percent_outside, percent_inside_i
for c_i in range(0, len(lon_lwind)):
lwind_contour_coords.append((lon_lwind[c_i], lat_lwind[c_i]))
lwind_contour = LineString(lwind_contour_coords)
ds_polar_north["low_wind_research"] = xr.DataArray(data=lwind_contour.wkt)
ds_polar_north["low_wind_research"] = xr.DataArray(data=lwind_contour.wkt,
attrs={"long_name": "WKT shape of area within which "
"low wind area is found"})
lon_hwind, lat_hwind = research_contour(msk_hwind, lon_algo, lat_algo)
hwind_contour_coords = []
for c_i in range(0, len(lon_lwind)):
hwind_contour_coords.append((lon_hwind[c_i], lat_hwind[c_i]))
hwind_contour = LineString(hwind_contour_coords)
ds_polar_north["high_wind_research"] = xr.DataArray(data=hwind_contour.wkt)
ds_polar_north["high_wind_research"] = xr.DataArray(data=hwind_contour.wkt,
attrs={"long_name": "WKT shape of area within which "
"high wind area is found"})
outfile = os.path.join(path, f"{os.path.basename(sat_file_sw).replace('.nc', '').replace('cm', 'ca')}_{sid}.nc")
logger.info(f"Writing netCDF product to {outfile}")
......
......@@ -35,7 +35,7 @@ def custom_meshgrid(theta, rad, data_shape):
coord_theta[i, j] = theta[i]
coord_rad[i, j] = rad[j]
print("Time taken ", time.time() - st)
#print("Time taken ", time.time() - st)
return coord_theta, coord_rad
......@@ -47,12 +47,12 @@ def to_polar(ds, lon_center, lat_center):
xloc, yloc = find_nearest_2d(ds.lon, ds.lat, lon_center, lat_center)
pix_dist = pixel_distance(ds, xloc, yloc)
print(pix_dist)
#print(pix_dist)
# Coordinates array
theta = np.linspace(0, 360, num=data_shape[0])
rads = np.linspace(0, max_radius, num=data_shape[1])
rads_meter = rads * pix_dist
print(rads_meter.min(), rads_meter.max())
#print(rads_meter.min(), rads_meter.max())
#coord_theta, coord_rad = np.meshgrid(theta, rads_meter)
......@@ -65,7 +65,7 @@ def to_polar(ds, lon_center, lat_center):
# Computing polar projection for all variables
vars_polar = {}
for var_name in ds.data_vars:
var_polar = cv2.linearPolar(ds[var_name].values, (int(yloc), int(xloc)), max_radius, cv2.WARP_FILL_OUTLIERS)
var_polar = cv2.linearPolar(ds[var_name].values, (int(yloc), int(xloc)), max_radius, cv2.INTER_NEAREST)
vars_polar[var_name] = (["theta", "rad"], var_polar)
vars_polar["lon"] = (["theta", "rad"], lons_polar)
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment