Good day. earlier this week i had an issue with processing daily precipitation and i think i solved it. below is the python code for extracting RAIN and RAINNC. The script extracts (RAINC and RAINNC and adds them to come with daily total)
import os
import xarray as xr
from datetime import datetime, timedelta
def increment_day(date_str):
# Function to increment the day by one
date_format = "%Y-%m-%d"
new_date = datetime.strptime(date_str, date_format) + timedelta(days=1)
return new_date.strftime(date_format)
def process_daily_rainfall(base_path, date_part, domain, precip_type):
# Format file names for the start of the current day and the next day
current_day_file = f'wrfout_{domain}_{date_part}_00:00:00'
next_day_file = f'wrfout_{domain}_{increment_day(date_part)}_00:00:00'
# Construct file paths
current_day_path = os.path.join(base_path, current_day_file)
next_day_path = os.path.join(base_path, next_day_file)
# Check if both required files are available
if not os.path.exists(current_day_path) or not os.path.exists(next_day_path):
print(f"Required files missing for date: {date_part}")
return
# Open datasets
with xr.open_dataset(current_day_path) as current_data, xr.open_dataset(next_day_path) as next_data:
try:
# Initialize variables for precipitation types
if precip_type == 'both':
daily_precip = (next_data['RAINC'] + next_data['RAINNC']) - (current_data['RAINC'] + current_data['RAINNC'])
elif precip_type == 'rainc':
daily_precip = next_data['RAINC'] - current_data['RAINC']
elif precip_type == 'rainnc':
daily_precip = next_data['RAINNC'] - current_data['RAINNC']
# Extract latitude and longitude (assuming they are the same for both datasets)
latitude = current_data['XLAT']
longitude = current_data['XLONG']
# Create a new dataset for the daily total
daily_total = xr.Dataset({
'TOTAL_PRECIP': daily_precip,
'XLAT': latitude,
'XLONG': longitude
})
# Save the dataset
output_file_path = os.path.join(base_path, f'wrfout_{domain}_{date_part}_daily_precip.nc')
daily_total.to_netcdf(output_file_path)
print(f'Daily total precipitation data for {date_part} saved to {output_file_path}')
except KeyError:
print("Error: Key variables RAINC or RAINNC missing in one of the files.")
def interactive_process():
base_path = '/data/run' # Update this path as needed for different directory structures
while True:
# Prompt for the date
date_part = input("Enter the date (YYYY-MM-DD) you want to process, or 'quit' to exit: ")
if date_part.lower() == 'quit':
print("Exiting the script.")
break
# Prompt for the domain (e.g., d01 or d02)
domain = input("Enter the domain (e.g., d01 or d02): ").strip()
# Prompt for precipitation type
precip_type = input("Select the precipitation type to process (RAINC, RAINNC, both): ").strip().lower()
# Validate precip_type input
if precip_type not in ['rainc', 'rainnc', 'both']:
print("Invalid precipitation type. Please enter RAINC, RAINNC, or both.")
continue
# Process the selected date and domain
process_daily_rainfall(base_path, date_part, domain, precip_type)
# Ask if the user wants to process another date
continue_choice = input("Do you want to process another date? (yes/no): ").strip().lower()
if continue_choice != 'yes':
print("Exiting the script.")
break
# Call the main interactive processing function
interactive_process()
Ideal Calculation for a Day (Day 2 Example)
- Required Data: To calculate the total rainfall for Day 2, you need:
- The cumulative rainfall at the end of Day 1 (this is the state at the beginning of Day 2, available in Day 2's file).
- The cumulative rainfall at the end of Day 2 (this is the state at the beginning of Day 3, available in Day 3's file).
- Calculation: The total rainfall for Day 2 is obtained by:
SUMMARY of calculation
Daily Rainfall on Day 2 = (Cumulative Rainfall at the end of Day 2) - (Cumulative Rainfall at the end of Day 1)
import os
import xarray as xr
from datetime import datetime, timedelta
def increment_day(date_str):
# Function to increment the day by one
date_format = "%Y-%m-%d"
new_date = datetime.strptime(date_str, date_format) + timedelta(days=1)
return new_date.strftime(date_format)
def process_daily_rainfall(base_path, date_part, domain, precip_type):
# Format file names for the start of the current day and the next day
current_day_file = f'wrfout_{domain}_{date_part}_00:00:00'
next_day_file = f'wrfout_{domain}_{increment_day(date_part)}_00:00:00'
# Construct file paths
current_day_path = os.path.join(base_path, current_day_file)
next_day_path = os.path.join(base_path, next_day_file)
# Check if both required files are available
if not os.path.exists(current_day_path) or not os.path.exists(next_day_path):
print(f"Required files missing for date: {date_part}")
return
# Open datasets
with xr.open_dataset(current_day_path) as current_data, xr.open_dataset(next_day_path) as next_data:
try:
# Initialize variables for precipitation types
if precip_type == 'both':
daily_precip = (next_data['RAINC'] + next_data['RAINNC']) - (current_data['RAINC'] + current_data['RAINNC'])
elif precip_type == 'rainc':
daily_precip = next_data['RAINC'] - current_data['RAINC']
elif precip_type == 'rainnc':
daily_precip = next_data['RAINNC'] - current_data['RAINNC']
# Extract latitude and longitude (assuming they are the same for both datasets)
latitude = current_data['XLAT']
longitude = current_data['XLONG']
# Create a new dataset for the daily total
daily_total = xr.Dataset({
'TOTAL_PRECIP': daily_precip,
'XLAT': latitude,
'XLONG': longitude
})
# Save the dataset
output_file_path = os.path.join(base_path, f'wrfout_{domain}_{date_part}_daily_precip.nc')
daily_total.to_netcdf(output_file_path)
print(f'Daily total precipitation data for {date_part} saved to {output_file_path}')
except KeyError:
print("Error: Key variables RAINC or RAINNC missing in one of the files.")
def interactive_process():
base_path = '/data/run' # Update this path as needed for different directory structures
while True:
# Prompt for the date
date_part = input("Enter the date (YYYY-MM-DD) you want to process, or 'quit' to exit: ")
if date_part.lower() == 'quit':
print("Exiting the script.")
break
# Prompt for the domain (e.g., d01 or d02)
domain = input("Enter the domain (e.g., d01 or d02): ").strip()
# Prompt for precipitation type
precip_type = input("Select the precipitation type to process (RAINC, RAINNC, both): ").strip().lower()
# Validate precip_type input
if precip_type not in ['rainc', 'rainnc', 'both']:
print("Invalid precipitation type. Please enter RAINC, RAINNC, or both.")
continue
# Process the selected date and domain
process_daily_rainfall(base_path, date_part, domain, precip_type)
# Ask if the user wants to process another date
continue_choice = input("Do you want to process another date? (yes/no): ").strip().lower()
if continue_choice != 'yes':
print("Exiting the script.")
break
# Call the main interactive processing function
interactive_process()