Opened 3 years ago
#309 new Defect
something wrong Matlab read the output "header" from FLEXPART10.4
Reported by: | zhangxun | Owned by: | |
---|---|---|---|
Priority: | major | Milestone: | FLEXPART 10 |
Component: | FP post-processing | Version: | FLEXPART 10.4 |
Keywords: | header | Cc: |
Description
After finished the backward run from FLEXPART 10.4,I used this code from:
Matlab functions to read FLEXPART gridded output and concentrations (S. Eckhardt, H. Sodemann, S. Raja) http://zardoz.nilu.no/~sabine/matlab4FLEXPART.tar.gz''',the function:
flex_header.m to read the header file. My code as followed:
###########################################
%% beginning clear all; clc; pathname = 'C:\Users\dell\Desktop\matlab4FLEXPART\'; % anonymous helper function to skip some bytes skip = @(fid) fread(fid,2,'int32'); nest=0; readp=1; % default: do not read release points calcarea=0; % datespath=pathname; header.path=pathname; fhd=fopen([header.path 'header']); %fprintf('\nreading header'); header.nest=0; %% % default settings header.decayconstant=0; % start reading if header was opened if fhd~=-1 rl=fread(fhd,1,'int32'); header.open=1; header.ibdate=fread(fhd,1,'int32'); header.ibtime=fread(fhd,1,'int32'); pivot_vec=(datevec(date)); header.simulationstart=datenum([num2str(header.ibdate) num2str(header.ibtime,'%06d')],'yyyymmddhhMMss'); header.ver=char(transpose(fread(fhd,13,'char'))); %% skip(fhd); header.loutstep=fread(fhd,1,'int32'); header.loutaver=fread(fhd,1,'int32'); header.loutsample=fread(fhd,1,'int32'); skip(fhd); header.outlon0=fread(fhd,1,'float'); header.outlat0=fread(fhd,1,'float'); header.numxgrid=fread(fhd,1,'int32'); header.numygrid=fread(fhd,1,'int32'); header.dxout=fread(fhd,1,'float'); header.dyout=fread(fhd,1,'float',8); header.numzgrid=fread(fhd,1,'int32'); header.outheight=fread(fhd,header.numzgrid,'float'); header.outheight=cat(1,[0],header.outheight); skip(fhd); header.jjjjmmdd=fread(fhd,1,'int32'); header.hhmmss=fread(fhd,1,'int32',8); header.nspec=fread(fhd,1,'int32')/3; end %%
###################################################
Here is my COMMAND file
*************************************************************************************************************** * * * Input file for the Lagrangian particle dispersion model FLEXPART * * Please select your options * * * *************************************************************************************************************** &COMMAND LDIRECT= -1, ! Simulation direction in time ; 1 (forward) or -1 (backward) IBDATE= 20120101, ! Start date of the simulation ; YYYYMMDD: YYYY=year, MM=month, DD=day IBTIME= 000000, ! Start time of the simulation ; HHMISS: HH=hours, MI=min, SS=sec; UTC IEDATE= 20120101, ! End date of the simulation ; same format as IBDATE IETIME= 230000, ! End time of the simulation ; same format as IBTIME LOUTSTEP= 3600, ! Interval of model output; average concentrations calculated every LOUTSTEP (s) LOUTAVER= 3600, ! Interval of output averaging (s) LOUTSAMPLE= 900, ! Interval of output sampling (s), higher stat. accuracy with shorter intervals ITSPLIT= 99999999, ! Interval of particle splitting (s) LSYNCTIME= 900, ! All processes are synchronized to this time interval (s) CTL= -5.0000000, ! CTL>1, ABL time step = (Lagrangian timescale (TL))/CTL, uses LSYNCTIME if CTL<0 IFINE= 4, ! Reduction for time step in vertical transport, used only if CTL>1 IOUT= 5, ! Output type: [1]mass 2]pptv 3]1&2 4]plume 5]1&4, +8 for NetCDF output IPOUT= 1, ! Particle position output: 0]no 1]every output 2]only at end 3]time averaged LSUBGRID= 1, ! Increase of ABL heights due to sub-grid scale orographic variations;[0]off 1]on LCONVECTION= 1, ! Switch for convection parameterization;0]off [1]on LAGESPECTRA= 0, ! Switch for calculation of age spectra (needs AGECLASSES);[0]off 1]on IPIN= 0, ! Warm start from particle dump (needs previous partposit_end file); [0]no 1]yes IOUTPUTFOREACHRELEASE= 1, ! Separate output fields for each location in the RELEASE file; [0]no 1]yes IFLUX= 0, ! Output of mass fluxes through output grid box boundaries MDOMAINFILL= 0, ! Switch for domain-filling, if limited-area particles generated at boundary IND_SOURCE= 2, ! Unit to be used at the source ; [1]mass 2]mass mixing ratio IND_RECEPTOR= 2, ! Unit to be used at the receptor; [1]mass 2]mass mixing ratio 3]wet depo. 4]dry depo. MQUASILAG= 0, ! Quasi-Lagrangian mode to track individual numbered particles NESTED_OUTPUT= 0, ! Output also for a nested domain LINIT_COND= 0, ! Output sensitivity to initial conditions (bkw mode only) [0]off 1]conc 2]mmr SURF_ONLY= 0, ! Output only for the lowest model layer, used w/ LINIT_COND=1 or 2 CBLFLAG= 0, ! Skewed, not Gaussian turbulence in the convective ABL, need large CTL and IFINE OHFIELDS_PATH= "../../flexin/", ! Default path for OH file /
#################################################
Here is my RELEASE file
*************************************************************************************************************** * * * * * * * Input file for the Lagrangian particle dispersion model FLEXPART * * Please select your options * * * * * * * *************************************************************************************************************** &RELEASES_CTRL NSPEC = 1, ! Total number of species SPECNUM_REL= 22, ! Species numbers in directory SPECIES / &RELEASE ! For each release IDATE1 = 20120101, ! Release start date, YYYYMMDD: YYYY=year, MM=month, DD=day ITIME1 = 000000, ! Release start time in UTC HHMISS: HH hours, MI=minutes, SS=seconds IDATE2 = 20120101, ! Release end date, same as IDATE1 ITIME2 = 230000, ! Release end time, same as ITIME1 LON1 = 126.542, ! Left longitude of release box -180 < LON1 <180 LON2 = 126.542, ! Right longitude of release box, same as LON1 LAT1 = 45.755, ! Lower latitude of release box, -90 < LAT1 < 90 LAT2 = 45.755, ! Upper latitude of release box same format as LAT1 Z1 = 0, ! Lower height of release box meters/hPa above reference level Z2 = 50.000, ! Upper height of release box meters/hPa above reference level ZKIND = 1, ! Reference level 1=above ground, 2=above sea level, 3 for pressure in hPa MASS = 1.0000E0, ! Total mass emitted, only relevant for fwd simulations PARTS = 10000, ! Total number of particles to be released COMMENT = "RELEASE 1", ! Comment, written in the outputfile / &RELEASE ! For each release IDATE1 = 20120101, ! Release start date, YYYYMMDD: YYYY=year, MM=month, DD=day ITIME1 = 000000, ! Release start time in UTC HHMISS: HH hours, MI=minutes, SS=seconds IDATE2 = 20120101, ! Release end date, same as IDATE1 ITIME2 = 230000, ! Release end time, same as ITIME1 LON1 = 126.561, ! Left longitude of release box -180 < LON1 <180 LON2 = 126.561, ! Right longitude of release box, same as LON1 LAT1 = 45.8167, ! Lower latitude of release box, -90 < LAT1 < 90 LAT2 = 45.8167, ! Upper latitude of release box same format as LAT1 Z1 = 0, ! Lower height of release box meters/hPa above reference level Z2 = 50.000, ! Upper height of release box meters/hPa above reference level ZKIND = 1, ! Reference level 1=above ground, 2=above sea level, 3 for pressure in hPa MASS = 1.0000E0, ! Total mass emitted, only relevant for fwd simulations PARTS = 10000, ! Total number of particles to be released COMMENT = "RELEASE 2", ! Comment, written in the outputfile / &RELEASE ! For each release IDATE1 = 20120101, ! Release start date, YYYYMMDD: YYYY=year, MM=month, DD=day ITIME1 = 000000, ! Release start time in UTC HHMISS: HH hours, MI=minutes, SS=seconds IDATE2 = 20120101, ! Release end date, same as IDATE1 ITIME2 = 230000, ! Release end time, same as ITIME1 LON1 = 126.979, ! Left longitude of release box -180 < LON1 <180 LON2 = 126.979, ! Right longitude of release box, same as LON1 LAT1 = 45.5422, ! Lower latitude of release box, -90 < LAT1 < 90 LAT2 = 45.5422, ! Upper latitude of release box same format as LAT1 Z1 = 0, ! Lower height of release box meters/hPa above reference level Z2 = 50.000, ! Upper height of release box meters/hPa above reference level ZKIND = 1, ! Reference level 1=above ground, 2=above sea level, 3 for pressure in hPa MASS = 1.0000E0, ! Total mass emitted, only relevant for fwd simulations PARTS = 10000, ! Total number of particles to be released COMMENT = "RELEASE 3", ! Comment, written in the outputfile /
#####################
Here is my header_txt file
# ibdate,ibtime, iedate, ietime, flexversion 20120101 0 20120101 230000 Version 10.4 (2019-11-12) # interval, averaging time, sampling time -3600 -3600 -900 # information on grid setup #outlon0,outlat0,numxgrid,numygrid,dxout,dyout 110.000000 40.0000000 200 200 0.100000001 0.100000001 # numzgrid, outheight(.) 5 10.0000000 100.000000 500.000000 1000.00000 2000.00000 # jjjjmmdd,ihmmss 20120101 230000 # information on species # 3*nspec,maxpointspec_act 3 3 # for nspec: # 1, WD_ # 1, DD_ # numzgrid,species 1 WD_CO 1 DD_CO 5 CO # information on model switches # method,lsubgrid,lconvection, ind_source,ind_receptor 0 1 1 2 2 # information on age class 1 999999999
###################################
Problem:
Actually, after Matlab read this header file,
header.loutstep should be "-3600", but the header strcture in matlab:
header.loutstep = 691155245.
header.loutaver should be "-3600", but the header strcture in matlab:
header.loutaver = 538976288
header.loutsample should be "-900", but the header strcture in matlab:
header.loutsample = 538976288
header.outlon0 shoule be "110.0000", but the header strcture is matlab:
header.outlon0 = 1.3563e-19.
What's more, header.outlat0, header.numxgrid, header.numygrid, header.dxout,
header.dyout ... the values are wrong!!!
header.jjjjmmdd, header.hhmmss, header.nspec, header.numpointspec in matlab
structure is NULL.
I wonder where is wrong when I used the matlab code read the header file values???
dates file