Skip to content

Commit

Permalink
Again, added features in bufr.py for carra2 to be able to read local …
Browse files Browse the repository at this point in the history
…observations
  • Loading branch information
PerDahlgren committed Jan 12, 2024
1 parent 213636c commit 3af9e20
Showing 1 changed file with 39 additions and 10 deletions.
49 changes: 39 additions & 10 deletions pysurfex/bufr.py
Original file line number Diff line number Diff line change
Expand Up @@ -37,6 +37,7 @@ def __init__(
latrange=None,
label="bufr",
use_first=False,
sigmao=None,
):
"""Initialize a bufr observation set.
Expand All @@ -49,6 +50,7 @@ def __init__(
latrange (list): Allowed range of latitides [min, max]
label (str): A label for the resulting observations set
use_first (bool): Use only the first valid observation for a point if more are found
sigmao (float, optional): Observation error relative to normal background error. Defaults to None.
Raises:
RuntimeError: ECCODES not found. Needed for bufr reading
Expand All @@ -66,6 +68,9 @@ def __init__(

# open bufr file
file_handler = open(bufrfile, mode="rb")
number_of_bytes = file_handler.seek(0, 2)
logging.info("File size: %s", number_of_bytes)
file_handler.seek(0)

# define the keys to be printed
keys = [
Expand All @@ -84,6 +89,7 @@ def __init__(
"blockNumber",
"stationOrSiteName",
]
processed_threshold = 0
nerror = {}
ntime = {}
nundef = {}
Expand Down Expand Up @@ -364,10 +370,28 @@ def __init__(
latrange[0] <= lat <= latrange[1]
and lonrange[0] <= lon <= lonrange[1]
):
obs_dtg = as_datetime_args(
year=year, month=month, day=day, hour=hour, minute=minute
)
if not np.isnan(value):
obs_dtg = None
try:
obs_dtg = as_datetime_args(
year=year,
month=month,
day=day,
hour=hour,
minute=minute,
)
except ValueError:
logging.warning(
"Bad observations time: year=%s month=%s day=%s hour=%s minute=%s Position is lon=%s, lat=%s",
year,
month,
day,
hour,
minute,
lon,
lat,
)
obs_dtg = None
if not np.isnan(value) and obs_dtg is not None:
if self.inside_window(obs_dtg, valid_dtg, valid_range):
logging.debug(
"Valid DTG for station %s %s %s %s %s %s %s %s",
Expand Down Expand Up @@ -407,18 +431,23 @@ def __init__(
ndomain.update({var: ndomain[var] + 1})

cnt += 1
try:
nbytes = file_handler.tell()
except ValueError:
nbytes = number_of_bytes

if (cnt % 1000) == 0:
print(".", end="")
sys.stdout.flush()
processed = int(round(float(nbytes) * 100.0 / float(number_of_bytes)))
if processed > processed_threshold and processed % 5 == 0:
processed_threshold = processed
logging.info("Read: %s%%", processed)

# delete handle
eccodes.codes_release(bufr)

logging.info("\nFound %s/%s", str(len(observations)), str(cnt))
logging.info("Found %s/%s", str(len(observations)), str(cnt))
logging.info("Not decoded: %s", str(not_decoded))
for var in variables:
logging.info("\nObservations for var=%s: %s", var, str(nobs[var]))
logging.info("Observations for var=%s: %s", var, str(nobs[var]))
logging.info(
"Observations removed because of domain check: %s", str(ndomain[var])
)
Expand All @@ -435,7 +464,7 @@ def __init__(
# close the file
file_handler.close()

ObservationSet.__init__(self, observations, label=label)
ObservationSet.__init__(self, observations, label=label, sigmao=sigmao)

@staticmethod
def td2rh(t_d, temp, kelvin=True):
Expand Down

0 comments on commit 3af9e20

Please sign in to comment.