Two-week HYSPLIT simulations of air parcel trajectories originating in Manitoba and transported across the globe, using ERA5 reanalysis data as meteorological input.
Tuesday, August 18, 2026
Video: Dust Storm
NOAA GOES-19 imagery showing a dust storm affecting the entire Winnipeg area on the afternoon of May 14, 2026.
Reference:
https://classic107.com/articles/rare-dust-storm-warning-issued-across-southern-manitoba
Thursday, August 6, 2026
Matlab: An example Map script
% Created by LI Xu % Version 1.0 % August 6, 2026 % If you have any question about this code, % please do not hesitate to contact me via E-mail: % jeremy456@163.com % Blog: % http://blog.sciencenet.cn/u/lixujeremy % http://lixuworld.blogspot.com/ clear; clc; close all; timebegin=tic; cur_data=date; cur_time=fix(clock); str1=sprintf('%s %.2d:%.2d:%.2d', cur_data, cur_time(4), cur_time(5), cur_time(6)); fprintf('Time Begin: '); fprintf(str1); fprintf('\n'); % Define the Google Street Map tile URL % name = 'google_streets'; name='streets'; % url = 'https://server.arcgisonline.com/ArcGIS/rest/services/World_Shaded_Relief/MapServer/tile/${z}/${y}/${x}'; % url = 'https://server.arcgisonline.com/ArcGIS/rest/services/World_Shaded_Relief/MapServer/tile/${z}/${y}/${x}'; % url = 'https://mt1.google.com/vt/lyrs=m&x=${x}&y=${y}&z=${z}&style=p.v:off|s.t:3|s.e:l|p.v:on'; % % url = 'https://mt1.google.com/vt/lyrs=m&x=${x}&y=${y}&z=${z}&style=p.v:off|s.t:3|s.e:l|p.v:on'; % url = 'https://mt1.google.com/vt/lyrs=m&x=${x}&y=${y}&z=${z}'; % url = 'https://server.arcgisonline.com/ArcGIS/rest/services/World_Physical_Map/MapServer/tile/${z}/${y}/${x}'; % attribution = 'Map data ©2026 Google'; % url = 'https://cartodb-basemaps-a.global.ssl.fastly.net/light_all/${z}/${x}/${y}.png'; % url = 'https://cartodb-basemaps-a.global.ssl.fastly.net/dark_all/${z}/${x}/${y}.png'; % name = 'carto_voyager'; % attribution = '© OpenStreetMap contributors, © CARTO'; % url = 'https://cartodb-basemaps-a.global.ssl.fastly.net/rastertiles/voyager/${z}/${x}/${y}.png'; % % % % Add it to your MATLAB basemap list % addCustomBasemap(name, url, 'Attribution', attribution); % url = 'https://server.arcgisonline.com/ArcGIS/rest/services/NatGeo_World_Map/MapServer/tile/${z}/${y}/${x}'; % url = 'https://server.arcgisonline.com/ArcGIS/rest/services/Specialty/NatGeo_World_Map/MapServer/tile/${z}/${y}/${x}'; % url = 'https://cartodb-basemaps-a.global.ssl.fastly.net/rastertiles/voyager/${z}/${x}/${y}.png'; attribution = '© OpenStreetMap contributors, © CARTO'; url = 'https://mt1.google.com/vt/lyrs=s&x=${x}&y=${y}&z=${z}'; addCustomBasemap(name, url, 'Attribution', attribution); % addCustomBasemap(name, url, 'Attribution', 'Sources: Esri, National Geographic, Garmin'); % url = 'https://stamen-tiles.a.ssl.fastly.net/toner/${z}/${x}/${y}.png'; % addCustomBasemap('stamen_toner', url, 'Attribution', 'Map tiles by Stamen Design, CC BY 3.0'); % Create geographic axes with OpenStreetMap basemap figure; gx = geoaxes('Basemap', name); set(gcf, 'OuterPosition', get(0, 'ScreenSize')); f = gcf; % Define custom latitude/longitude range (Winnipeg-centered) % latRange = [49.76, 50.01]; % lonRange = [-97.45, -96.9]; latRange = [49.79+0.005+0.009, 49.82-0.0068]; lonRange = [-97.15+0.008, -97.13+0.004]; % Calculate aspect ratio for exact limits without expansion mean_lat = mean(latRange); cos_factor = cosd(mean_lat); ar = diff(latRange) / (diff(lonRange) * cos_factor); % Required height/width for axes % Set up resize function to maintain aspect and exact limits gx.Units = 'pixels'; f.Units = 'pixels'; f.SizeChangedFcn = @(~, ~) myresizefcn(gx, f, ar, latRange, lonRange); myresizefcn(gx, f, ar, latRange, lonRange); % Rotate the whole map figure visually view(gx, 2); % rotate3d(f, 'on'); % % rotate3d(gx, 'on'); % gx.View = [30 90]; gx.Grid = 'off'; gx.TickDir = 'none'; gx.LatitudeLabel.String = ''; gx.LongitudeLabel.String = ''; gx.FontSize = 0.01; % dat=cell2mat(dat); % lats=dat(:, 1); % lons=dat(:, 2); % nums=dat(:, 3); % hold on; % for i = 1:length(lats) % % Draw Circle % plot(gx, lats(i), lons(i), 'bo', 'MarkerSize', 18, 'MarkerFaceColor', 'y'); % % Draw Number % text(gx, lats(i), lons(i), num2str(nums(i)), ... % 'HorizontalAlignment', 'center', 'VerticalAlignment', 'middle', 'Color', 'blue', 'FontWeight','bold', 'FontSize', 15); % end % 49.81155561273022, -97.13142673818446 plot(gx, 49.81155561273022, -97.13142673818446, 'pentagram', 'MarkerSize', 20, 'MarkerFaceColor', 'blue'); hold off; gx.Grid = 'off'; gx.TickDir = 'none'; gx.LatitudeLabel.String = ''; gx.LongitudeLabel.String = ''; gx.FontSize = 0.01; % dim = [0.1 0.83 0.8 0.1]; % Position at the top center % str = 'Campus Water Fountain Map'; % % hAnn = annotation('textbox', dim, 'String', str, ... % 'EdgeColor', 'none', ... % No border % 'HorizontalAlignment', 'center', ... % 'VerticalAlignment', 'middle', ... % 'FontSize', 25, ... % 'FontWeight', 'bold', ... % 'Color', 'white', ... % 'Interpreter', 'none'); pause(3); otname='kunakmap_vector'; export_fig(otname, '-png', '-painters', '-r300'); close all; fprintf('Time Begin: '); fprintf(str1); fprintf('\n'); cur_data=date; cur_time=fix(clock); str2=sprintf('%s %.2d:%.2d:%.2d', cur_data, cur_time(4), cur_time(5), cur_time(6)); fprintf('Time End: '); disp(str2); timespan=toc(timebegin); fprintf('Time Span: %.4f s\n', timespan); disp('************************************************************'); function myresizefcn(gx, f, ar, latlim, lonlim) if (f.Position(4) / f.Position(3)) < ar % Figure is wide: constrain by height gx.Position(4) = f.Position(4) * 0.8; gx.Position(3) = gx.Position(4) / ar; % Center gx.Position(1) = (f.Position(3) - gx.Position(3)) / 2 + f.Position(1); gx.Position(2) = (f.Position(4) - gx.Position(4)) / 2 + f.Position(2); else % Figure is tall: constrain by width gx.Position(3) = f.Position(3) * 0.8; gx.Position(4) = gx.Position(3) * ar; % Center gx.Position(1) = (f.Position(3) - gx.Position(3)) / 2 + f.Position(1); gx.Position(2) = (f.Position(4) - gx.Position(4)) / 2 + f.Position(2); end % Re-apply limits to maintain exact range geolimits(gx, latlim, lonlim); end
Tuesday, May 12, 2026
Matlab: Table on the PPT created via Matlab
% Created by LI Xu % Version 1.0 % August 28, 2025 % If you have any question about this code, % please do not hesitate to contact me via E-mail: % jeremy456@163.com % Blog: % http://blog.sciencenet.cn/u/lixujeremy % http://lixuworld.blogspot.com/ clear; clc; close all; timebegin=tic; cur_data=date; cur_time=fix(clock); str1=sprintf('%s %.2d:%.2d:%.2d', cur_data, cur_time(4), cur_time(5), cur_time(6)); fprintf('Time Begin: '); fprintf(str1); fprintf('\n'); import mlreportgen.ppt.* % 1. Create a presentation ppt = Presentation('AQ.pptx'); % 2. Create a slide (e.g., 'Blank' or 'Title and Content') slide = add(ppt, 'Title and Content'); replace(slide, 'Title', 'Measured Data Table'); % 3. Prepare your data (Cell array: Rows x Cols) data = {'Variable', 'Campus', 'Winnipeg'; 'Temperature (℃) ', '40', '15'; ... 'Relative Humidity (%)', '65', '70'; ... 'Wind speed (m/s)', '65', '70'; ... 'PM2.5 (μg/m³)', '65', '70'; ... 'PM10 (μg/m³)', '65', '70'; ... 'NO₂ (ppb)', '65', '70'; ... 'O₃ (ppb)', '65', '70'; ... 'CO₂ (ppm)', '65', '70'}; % 4. Create and add the table % 'Table' automatically creates rows/columns based on the data tableObj = Table(data); firstRow = tableObj.Children(1); firstRow.Style = {BackgroundColor('#296080'), FontColor('white'), Bold()}; tableObj.Style = {HAlign('center'), VAlign('middle'), Border('solid'), RowHeight('0.3in')}; % tableObj.ColWidths = {'1.5in', '1.0in', '1.0in'}; % Column 3: 1.0 inches (Winnipeg) spec1 = ColSpec('2.5 in'); spec2 = ColSpec('1.5 in'); spec3 = ColSpec('1.5 in'); % 2. Assign the array of specs to the table tableObj.ColSpecs = [spec1, spec2, spec3]; % 1. Target the specific cell (Row 2, Column 2) using entry(row, col) targetCell = entry(tableObj, 2, 2); % 2. Apply the background color targetCell.Style = {BackgroundColor('red')}; % 1. Target the specific cell (e.g., Row 2, Column 2) targetCell = entry(tableObj, 2, 2); % 2. Apply the font color (e.g., red) % You can use names like 'red', 'white', or hex codes like '#FF0000' targetCell.Style = [targetCell.Style, {FontColor('white')}]; % Add the table to the content placeholder replace(slide, 'Content', tableObj); % 5. Close and View close(ppt); rptview(ppt); fprintf('Time Begin: '); fprintf(str1); fprintf('\n'); cur_data=date; cur_time=fix(clock); str2=sprintf('%s %.2d:%.2d:%.2d', cur_data, cur_time(4), cur_time(5), cur_time(6)); fprintf('Time End: '); disp(str2); timespan=toc(timebegin); fprintf('Time Span: %.4f s\n', timespan); disp('*************************************************************************');
Thursday, April 9, 2026
Data: Manitoba Ag-Weather Network Weather Stations
Manitoba Agriculture operates a network of over 100 weather stations across Manitoba's agricultural region which provide hourly updated air temperature, relative humidity, precipitation, wind speed and direction, soil temperature and soil moisture. The WeatherStations feature class displays the location of the weather stations and its table contains a link to current weather for each station. list
How to download
How to use them:
- Select a Station Name from the dropdown (they match the ones in your CSV).
- Pick a Date (or planting/end dates for seasonal).
- Submit to view the table for that station/date.
- Data goes back several years (stations have varying start dates, many from ~2011 onward, with 15-minute raw intervals rolled up to hourly/daily). Note: Raw data can have occasional sensor issues — use with caution.
Other useful pages:
- Current + summary pages: https://www.gov.mb.ca/agriculture/weather/weather-conditions-and-reports.html (links to the reports above + soil temps).
- Manitoba Ag Weather Program overview: https://www.gov.mb.ca/agriculture/weather/manitoba-ag-weather.html (confirms data is free and archived).
- Interactive current map: From the original dashboard or https://www.gov.mb.ca/agriculture/weather/current-ag-weather-conditions.html.
Code.
Thursday, April 2, 2026
Wednesday, March 18, 2026
Hysplit: To test Dry Deposition of Gases in HYSPLIT using a built-in example
- Prepare the Basic Run Setup (if not already done)
- Open the main Setup Run window.
- Set Total run time (hrs): e.g., 12 or 25 (tutorial uses 25 h for noticeable deposition).
- Set Direction: Fwd (forward).
- Add meteorology: Click Add Meteorology Files and select your file (e.g., oct1618.BIN).
- Click Setup starting locations → define 1 source (e.g., lat 40.0, lon -90.0, height 50 m, emission rate suitable for test).
- Set Pollutant & Deposition (Gas Dry Deposition)
- In main setup → click Pollutant, Concentration Grid, and Deposition setup.
- Under Pollutant: Num=1, select Specie 1.
- Under Deposition: Num=1, select Specie 1 → this opens the Deposition Definition window.
- In Deposition window:
- Particle or Gas: Gas.
- Dry Deposition: Yes.
- Wet Deposition: No.
- Preconfigure: Select I131g (built-in gaseous I-131 example; auto-loads realistic parameters like Vd ~0.006 m/s, molecular weight, Henry's law constants).
- Or manually: Set Vel(m/s) to 0.01 (1 cm/s, tutorial value for strong test).
- Leave other fields default/zero unless needed.
- Click green Save → close window.
- Add Ground Deposition Level to Concentration Grid
- In the same Pollutant, Concentration Grid, and Deposition setup window:
- Under Grids: Num=1 → click radio for Grid 1 (or double-click/label) → opens Definition of Concentration Grid 1.
- In Grid Definition window:
- Num of vertical levels: Change from 1 to 2.
- Height of levels(M Agl): Change to 0 100 (space-separated; 0 first for deposition accumulation, 100 for air concentration layer).
- Optional: Update Center of Lat and Lon to match your source (e.g., 40.0 -90.0).
- Keep spacing/span reasonable (e.g., 0.05 deg, 30 deg span).
- Output grid file name: Keep cdump or change.
- Click green Save → close.
- Run the Model
- Back in main Concentration Setup window → Save.
- Click Run Model (or equivalent button to execute HYSPLIT).
- Wait for completion → open the MESSAGE file in your working directory (text file named MESSAGE):
- Note the final total mass (near end, often in a NOTICE or summary line).
- Compare to initial emitted mass (early in file or from your emission setup: rate × duration).
- Look for vertical mass distribution tables (e.g., %Mass by height layer) — reduced low-level % indicates dry removal to ground.
- Compute % mass lost = (initial - final) / initial × 100. Tutorial expects ~10-14% for Vd=0.01 m/s over longer runs; shorter runs or lower Vd (like 0.006) show less.
- This confirms dry deposition occurred (mass removed from atmosphere).
- Display the Results (Air Concentration + Deposition Footprint)
- Go to Display → Concentration Contours (opens your Concentration Display window).
- Key settings:
- Input File: cdump (or your output name).
- Select Pollutant: All or your species (e.g., I131 / II131).
- Vertical Display: Show Each Level.
- From Bottom Level: 0 (deposition/surface).
- Through Top Level: 100 (air layer).
- Deposition Multiplier: Select Total (shows cumulative/total deposited mass over entire run; critical for seeing the full footprint instead of incremental amounts).
- Keep numeric field 1.0 (unless scaling units).
- Concentration Multiplier: 1.0.
- Uncheck Exposure if checked.
- Contour drawing options: Dyn-Exp (dynamic exponential; best for plumes).
- Label Source Rings: On (shows distance rings around source; set Number=4, Dist=100 km, Center to your source lat/lon).
- Output File: e.g., dry_depo_plot (unique name).
- View On: Checked.
- Click green Execute Display.
- What you'll see (sequence through frames):
- Air concentration frames (at 100 m level): Reduced plume due to deposition loss.
- Final/cumulative frame (level 0 with Total): Colored deposition footprint on ground (accumulated mass deposited as plume passed).
- If blank: Confirm level 0 in grid, Vd > 0, re-run if needed.
Subscribe to:
Posts (Atom)