@@ -106,6 +106,33 @@ def get_fraction_per_soil_type(soil_id):
106106 return soil_fractions_target
107107
108108
109+ def calculate_soil_fraction_optimized (target_grid , soil_types_raw , nearest_target_cell_to_raw_cells , ncpu = 13 ):
110+ """
111+ target_grid: target ICON grid
112+ soil_types_raw: landuse class for each cell from the HWSD dataset (LU variable)
113+ nearest_target_cell_to_raw_cell: indices of the cell from the target ICON grid which is nearest to each cell of the raw grid (from HWSD dataset)
114+ """
115+ ncells_target = target_grid .lons .size
116+ nsoil_types = 13
117+ nthreads = min (nsoil_types , ncpu )
118+
119+ soil_fractions_target = np .zeros ((ncells_target , nsoil_types ))
120+
121+ n_nearest_raw_cells = np .bincount (nearest_target_cell_to_raw_cells .ravel (),
122+ minlength = ncells_target )
123+
124+ valid_mask = (soil_types_raw >= 1 ) & (soil_types_raw <= nsoil_types )
125+ counts = np .zeros ((ncells_target , nsoil_types ), dtype = int )
126+ np .add .at (counts , (nearest_target_cell_to_raw_cells [valid_mask ], soil_types_raw [valid_mask ].astype (int )- 1 ), 1 )
127+
128+ np .divide (counts ,
129+ n_nearest_raw_cells [:, np .newaxis ],
130+ out = soil_fractions_target ,
131+ where = n_nearest_raw_cells [:, np .newaxis ] != 0 )
132+
133+ return soil_fractions_target
134+
135+
109136# --------------------------------------------------------------------------
110137# --------------------------------------------------------------------------
111138# initialize logger
0 commit comments