From f632ea35aa16e424b0cfef0fc240d3df31e1719a Mon Sep 17 00:00:00 2001 From: Felix Kleinert <f.kleinert@fz-juelich.de> Date: Wed, 31 Mar 2021 17:18:52 +0200 Subject: [PATCH] prepare sector assingmnet --- mlair/data_handler/data_handler_wrf_chem.py | 72 ++++++++++++++++++--- 1 file changed, 62 insertions(+), 10 deletions(-) diff --git a/mlair/data_handler/data_handler_wrf_chem.py b/mlair/data_handler/data_handler_wrf_chem.py index c7313ba6..176fdd8f 100644 --- a/mlair/data_handler/data_handler_wrf_chem.py +++ b/mlair/data_handler/data_handler_wrf_chem.py @@ -548,16 +548,17 @@ class DataHandlerSingleGridColumn(DataHandlerSingleStation): self.make_history_window(self.target_dim, self.window_history_size, self.time_dim) if self.var_logical_z_coord_selector is not None: self.history = self.history.sel({self._logical_z_coord_name: self.var_logical_z_coord_selector}) - self.history = self.modify_history() self.make_labels(self.target_dim, self.target_var, self.time_dim, self.window_lead_time) self.make_observation(self.target_dim, self.target_var, self.time_dim) if self.targetvar_logical_z_coord_selector is not None: self.label = self.label.sel({self._logical_z_coord_name: self.targetvar_logical_z_coord_selector}) self.observation = self.observation.sel( {self._logical_z_coord_name: self.targetvar_logical_z_coord_selector}) + + self.remove_nan(self.time_dim) + self.history = self.modify_history() self.label = self.modify_label() self.observation = self.modify_observation() - self.remove_nan(self.time_dim) # self.input_data = self.input_data.compute() # self.label = self.label.compute() @@ -686,11 +687,11 @@ class DataHandlerSectorGrid(DataHandlerSingleGridColumn): trafo_opts = self._transformation[pos_of_trafo] return trafo_opts - def apply_transformation_on_wind_sector_edges(self, ws_edges): - ws_trafo_opts = self.get_transformation_opts() - ws_trafo, _ = self.transform(ws_edges, dim=self.aggregation_dim, inverse=False, - opts=ws_trafo_opts, transformation_dim=self.target_dim) - return ws_trafo + def apply_transformation_on_data(self, data): + trafo_opts = self.get_transformation_opts() + trafo, _ = self.transform(data, dim=self.aggregation_dim, inverse=False, + opts=trafo_opts, transformation_dim=self.target_dim) + return trafo def modify_history(self): if self.transformation_is_applied: @@ -699,12 +700,63 @@ class DataHandlerSectorGrid(DataHandlerSingleGridColumn): sector_allocation = self.windsector.get_sect_of_value(value=wind_dir_of_interest, external_edges=ws_edges) existing_sectors = np.unique(sector_allocation.data) with self.loader as loader: - pass - #circular_data = loader.data[self.variables].where(loader.geo_infos.dist.squeeze() <= self.radius) + # setup sector history + sector_history = xr.ones_like(self.history) + sector_history_var_names = [f"{var}_sect" for var in sector_history.coords[self.target_dim].values] + # sector_history = sector_history.assign_coords({self.target_dim: sector_history_var_names}) + + grid_data = self.preselect_and_transform_neighbouring_data_based_on_radius(loader) + + for sect in existing_sectors: + # select data in wind sector + sec_data = self.get_section_data_from_circle(grid_data, loader, sect) + sec_data = self.apply_aggregation_method_on_sector_data(sec_data, loader) + + # loader.data.T2.where(loader.geo_infos.dist.sel({'points': 0}).drop('points') <= self.radius).where( + # self.windsector.is_in_sector(sect, loader.geo_infos.bearing)) + # + # loader.data[self.variables].sel( + # {self.time_dim: self.history.coords[self.time_dim].values}).where( + # loader.geo_infos.dist.sel({'points': 0}).drop('points') <= self.radius).where( + # self.windsector.is_in_sector(sect, + # loader.geo_infos.bearing.drop('points').squeeze())) return self.history + def get_section_data_from_circle(self, grid_data, loader, sect): + sec_data = grid_data.where( + self.windsector.is_in_sector(sect, loader.geo_infos.bearing.drop('points').squeeze())) + return sec_data + + def preselect_and_transform_neighbouring_data_based_on_radius(self, loader): + """ + Select neighbouring grid boxes which have a maximal distance of pre-selected radius from full model field + :param loader: + :type loader: + :return: + :rtype: + """ + # get data from loader + grid_data = loader.data[self.variables] + # select correct time steps and vertical layers + grid_data = grid_data.sel({ + # self.time_dim: self.history.coords[self.time_dim].values, # get data which is used in history + self.time_dim: slice(self.start, self.end), # get data which is used in history + self._logical_z_coord_name: self.var_logical_z_coord_selector, # get only used vertical layers + }) + # select grid boxes which are closer than given radius + grid_data = grid_data.where( + loader.geo_infos.dist.sel({"points": 0}).drop("points") <= self.radius) + # apply transformation on variables + grid_data = self.apply_transformation_on_data(grid_data.to_array(self.target_dim)).to_dataset(self.target_dim) + + return grid_data + + def apply_aggregation_method_on_sector_data(self, data, loader): + data = data.mean(dim=(loader.logical_x_coord_name, loader.logical_y_coord_name)) + return data + def compute_wind_dir_of_interest(self): wind_dir_of_intrest = self.history.sel({self.target_dim: self.wind_dir_name, self.window_dim: 0}) return wind_dir_of_intrest @@ -712,7 +764,7 @@ class DataHandlerSectorGrid(DataHandlerSingleGridColumn): @TimeTrackingWrapper def get_applied_transdormation_on_wind_sector_edges(self): ws_edges = self._get_left_and_right_wind_sector_edges(return_as='xr.da', dim=self.wind_sector_edge_dim_name) - ws_edges = self.apply_transformation_on_wind_sector_edges(ws_edges) + ws_edges = self.apply_transformation_on_data(ws_edges) return ws_edges # def set_inputs_and_targets(self): -- GitLab