1 | ;+--------------------------------------------------------------------------- |
---|
2 | PRO read_cdfobs, Files, NumObs=NumObs, $ |
---|
3 | Latitudes=Latitudes, Longitudes=Longitudes, Depths=Depths, $ |
---|
4 | Obs=Observations, ModelVals=ModelVals, qcs=QCs, $ |
---|
5 | Ob2=Observations2, ModelVal2=ModelVals2, qc2=QCs2, $ |
---|
6 | Ob3=Observations3, $ |
---|
7 | Dates=Dates, rmdi=rmdi, iobs=iobs, jobs=jobs, kobs=kobs, $ |
---|
8 | Salinity=Salinity, nodates=nodates, types=types, $ |
---|
9 | filetype=ObsType, quiet=quiet, MDT=MDT, $ |
---|
10 | ProfileNum=ProfileNum, error=error, $ |
---|
11 | notfussy=notfussy,VarName=VarName |
---|
12 | ;+-------------------------------------------------------------------------- |
---|
13 | ; Read in observation and feedback files |
---|
14 | ; detects filetype and calls the appropriate reading routine |
---|
15 | ; |
---|
16 | ; Author: D. J. Lea Feb 2008 |
---|
17 | ;+-------------------------------------------------------------------------- |
---|
18 | |
---|
19 | ; Declare error label. |
---|
20 | ON_IOERROR, IOERROR |
---|
21 | |
---|
22 | error=0 |
---|
23 | title='' |
---|
24 | |
---|
25 | ; Set types to undefined |
---|
26 | types=-1 & tempvar=size(temporary(types)) |
---|
27 | |
---|
28 | ; Read in netcdf data file and observation operator |
---|
29 | ;2. Work out which type of data is in the files by looking at the first one. |
---|
30 | ncid = ncdf_open(Files(0), /nowrite) |
---|
31 | if ncid lt 0 then message, 'Error opening file '+File |
---|
32 | result = NCDF_ATTINQ( ncid, 'title', /global) |
---|
33 | if result.datatype ne 'UNKNOWN' then ncdf_attget, ncid, 'title', Title, /global |
---|
34 | |
---|
35 | if string(Title) eq "NEMO observation operator output" then ObsType='feedback' $ |
---|
36 | else $ |
---|
37 | if string(Title) eq "Forecast class 4 file" then ObsType='ForecastClass4' $ |
---|
38 | else $ |
---|
39 | ObsType = 'none' |
---|
40 | |
---|
41 | if ObsType eq 'none' then begin |
---|
42 | varid = ncdf_varid(ncid, 'POTM_CORRECTED') |
---|
43 | if varid ge 0 then ObsType = 'Prof' |
---|
44 | |
---|
45 | varid = ncdf_varid(ncid, 'SLA') |
---|
46 | if varid ge 0 then ObsType = 'SSH' |
---|
47 | |
---|
48 | varid = ncdf_varid(ncid, 'SST') |
---|
49 | varid2 = ncdf_varid(ncid, 'sea_surface_temperature') |
---|
50 | varid=max([varid,varid2]) |
---|
51 | if varid ge 0 then ObsType = 'SST' |
---|
52 | |
---|
53 | varid = ncdf_varid(ncid, 'SEAICE') |
---|
54 | varid2 = ncdf_varid(ncid, 'sea_ice_concentration') |
---|
55 | varid=max([varid,varid2]) |
---|
56 | if varid ge 0 then ObsType = 'SEAICE' |
---|
57 | |
---|
58 | varid = ncdf_varid(ncid, 'LOGCHL') |
---|
59 | varid2 = ncdf_varid(ncid, 'LogChl') |
---|
60 | varid3 = ncdf_varid(ncid, 'CHL1_mean') |
---|
61 | varid=max([varid,varid2,varid3]) |
---|
62 | if varid ge 0 then ObsType = 'LOGCHL' |
---|
63 | |
---|
64 | varid = ncdf_varid(ncid, 'PCO2') |
---|
65 | if varid ge 0 then ObsType = 'PCO2' |
---|
66 | |
---|
67 | varid = ncdf_varid(ncid, 'UCRT') |
---|
68 | varid2 = ncdf_varid(ncid, 'VCRT') |
---|
69 | varid=max([varid,varid2]) |
---|
70 | if varid ge 0 then ObsType = 'CRT' |
---|
71 | endif |
---|
72 | |
---|
73 | ncdf_close, ncid |
---|
74 | |
---|
75 | if not keyword_set(quiet) then print, 'Reading in data as type ',ObsType |
---|
76 | |
---|
77 | if ObsType eq 'feedback' then begin |
---|
78 | |
---|
79 | if n_elements(DepRange) eq 0 then DepRange=[0,5000] |
---|
80 | |
---|
81 | read_feedback, Files, DepRange=DepRange, VarName=VarName, NumData=numobs, $ |
---|
82 | OutLats=Latitudes, OutLons=Longitudes, $ |
---|
83 | OutInstruments=Instruments, OutPlatform=Platform, $ |
---|
84 | OutDeps=Depths, $ |
---|
85 | OutObs1=Observations, OutObs2=Observations2, $ |
---|
86 | OutMod1=ModelVals, OutMod2=ModelVals2, $ |
---|
87 | OutQC1=QCs, OutQC2=QCs2, $ |
---|
88 | MDT=MDT, OutDates=Dates, $ |
---|
89 | rmdi=rmdi, quiet=quiet, $ |
---|
90 | OutProfileNum=ProfileNum |
---|
91 | |
---|
92 | ; info, instruments |
---|
93 | ; info, platform |
---|
94 | |
---|
95 | types=Platform |
---|
96 | if n_elements(instruments) gt 0 then types=Platform+' '+Instruments |
---|
97 | |
---|
98 | print, 'Varname: ',VarName |
---|
99 | endif else if ObsType eq 'Prof' then begin |
---|
100 | |
---|
101 | if keyword_set(PlotModel) then begin |
---|
102 | OutTMod = 1 |
---|
103 | OutSMod = 1 |
---|
104 | endif |
---|
105 | |
---|
106 | if n_elements(DepRange) eq 0 then DepRange=[0,5000] |
---|
107 | |
---|
108 | read_enact, Files, DepRange=DepRange, NumData=NumData, $ |
---|
109 | OutLats=Latitudes, OutLons=Longitudes, $ |
---|
110 | Instruments=Instruments, Platform=Platform, $ |
---|
111 | OutDeps=Depths, OutSObs=OutSObs, OutSQC=OutSQC, $ |
---|
112 | OutTObs=OutTObs, OutTQC=OutTQC, OutDates=Dates, $ |
---|
113 | OutTMod=OutTMod, OutSMod=OutSMod, rmdi=rmdi, quiet=quiet, $ |
---|
114 | iobs=iobs, jobs=jobs, kobs=kobs, $ |
---|
115 | ProfileNum=ProfileNum |
---|
116 | |
---|
117 | types=Instruments+' '+Platform |
---|
118 | |
---|
119 | ; if salinity keyword is set then put salinity values first |
---|
120 | if keyword_set(Salinity) then begin |
---|
121 | Observations = OutSObs |
---|
122 | ModelVals = OutSMod |
---|
123 | QCs = OutSQC |
---|
124 | Observations2 = OutTObs |
---|
125 | ModelVals2 = OutTMod |
---|
126 | QCs2 = OutTQC |
---|
127 | endif else begin |
---|
128 | Observations = OutTObs |
---|
129 | ModelVals = OutTMod |
---|
130 | QCs = OutTQC |
---|
131 | Observations2 = OutSObs |
---|
132 | ModelVals2 = OutSMod |
---|
133 | QCs2 = OutSQC |
---|
134 | endelse |
---|
135 | |
---|
136 | numobs=numdata |
---|
137 | numobs=n_elements(latitudes) |
---|
138 | |
---|
139 | endif else if ObsType eq 'SSH' then begin |
---|
140 | |
---|
141 | read_sla, Files, NumObs=NumObs, $ |
---|
142 | Latitudes=Latitudes, Longitudes=Longitudes, $ |
---|
143 | ObsSLA=Observations, ModSLA=ModelVals, SLAQC=QCs, $ |
---|
144 | MDT=MDT, Satellites=Satellites, types=types, Dates=Dates, rmdi=rmdi, $ |
---|
145 | iobs=iobs, jobs=jobs, mstp=mstp, quiet=quiet |
---|
146 | Depths = fltarr(NumObs) |
---|
147 | |
---|
148 | endif else if ObsType eq 'SST' then begin |
---|
149 | |
---|
150 | read_sst, Files, NumObs=NumObs, $ |
---|
151 | Latitudes=Latitudes, Longitudes=Longitudes, $ |
---|
152 | ObsSST=Observations, ModSST=ModelVals, SSTQC=QCs, $ |
---|
153 | Dates=Dates, rmdi=rmdi, types=types, iobs=iobs, jobs=jobs, $ |
---|
154 | quiet=quiet |
---|
155 | Depths = fltarr(NumObs) |
---|
156 | |
---|
157 | endif else if ObsType eq 'SEAICE' then begin |
---|
158 | |
---|
159 | read_seaice, Files, NumObs=NumObs, $ |
---|
160 | Latitudes=Latitudes, Longitudes=Longitudes, $ |
---|
161 | Obs=Observations, Modarr=ModelVals, QC=QCs, $ |
---|
162 | Dates=Dates, rmdi=rmdi, iobs=iobs, jobs=jobs, nodates=nodates, $ |
---|
163 | quiet=quiet |
---|
164 | Depths = fltarr(NumObs) |
---|
165 | |
---|
166 | endif else if ObsType eq 'LOGCHL' then begin |
---|
167 | |
---|
168 | read_chl, Files, NumObs=NumObs, $ |
---|
169 | Latitudes=Latitudes, Longitudes=Longitudes, $ |
---|
170 | Obs=Observations, Modarr=ModelVals, QC=QCs, $ |
---|
171 | Dates=Dates, rmdi=rmdi, iobs=iobs, jobs=jobs, nodates=nodates, $ |
---|
172 | quiet=quiet |
---|
173 | Depths = fltarr(NumObs) |
---|
174 | |
---|
175 | endif else if ObsType eq 'PCO2' then begin |
---|
176 | |
---|
177 | read_pco2, Files, NumObs=NumObs, $ |
---|
178 | Latitudes=Latitudes, Longitudes=Longitudes, $ |
---|
179 | Obs=Observations, Modarr=ModelVals, QC=QCs, $ |
---|
180 | Dates=Dates, rmdi=rmdi, iobs=iobs, jobs=jobs, nodates=nodates, $ |
---|
181 | quiet=quiet |
---|
182 | Depths = fltarr(NumObs) |
---|
183 | |
---|
184 | endif else if ObsType eq 'CRT' then begin |
---|
185 | |
---|
186 | read_crt, Files, NumObs=NumObs, $ |
---|
187 | OutLats=Latitudes, OutLons=Longitudes, $ |
---|
188 | OutTypes=Types, OutDates=Dates, $ |
---|
189 | OutUObs=Observations, OutVObs=Observations2, $ |
---|
190 | OutUMod=ModelVals, OutVMod=ModelVals2, Quiet=Quiet, $ |
---|
191 | FloatNum=FloatNum, QC=QCs, OutSpeed=Observations3 |
---|
192 | rmdi=0.0 |
---|
193 | QCs2=QCs |
---|
194 | Depths = fltarr(NumObs) |
---|
195 | |
---|
196 | endif else if ObsType eq 'ForecastClass4' then begin |
---|
197 | |
---|
198 | print, 'reading forecast class 4 files' |
---|
199 | read_forc_class4, Files, NumObs=NumObs, $ |
---|
200 | Latitudes=Latitudes, Longitudes=Longitudes, $ |
---|
201 | Depths=Depths, $ |
---|
202 | Types=types, Dates=Dates, $ |
---|
203 | Obsarr=Observations, Obs2arr=Observations2, $ |
---|
204 | Modarr=ModelVals, Mod2arr=ModelVals2, $ |
---|
205 | QCs1=QCs, QCs2=QCs2, rmdi=rmdi, $ |
---|
206 | notfussy=notfussy |
---|
207 | |
---|
208 | endif else message, 'Error: ObsType is not set correctly' |
---|
209 | |
---|
210 | if (n_elements(types) eq 0) then begin |
---|
211 | types=replicate(rmdi,NumObs) |
---|
212 | endif |
---|
213 | |
---|
214 | goto, NOERROR |
---|
215 | |
---|
216 | IOERROR: error=1 |
---|
217 | print,'read_cdfobs: an error occurred trying read files: ',files |
---|
218 | |
---|
219 | NOERROR: |
---|
220 | |
---|
221 | END |
---|