Changeset 175 for trunk/src/sst_correction_ncdf.pro
- Timestamp:
- 03/22/12 14:40:17 (12 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/src/sst_correction_ncdf.pro
r174 r175 1 1 ;+ 2 ;3 ; .. _sst_correction_ncdf.pro:4 2 ; 5 3 ; ======================= … … 7 5 ; ======================= 8 6 ; 7 ; .. function:: sst_correction_ncdf(yyyymmddb,yyyymmdde) 8 ; 9 ; DESCRIPTION 10 ; =========== 11 ; 9 12 ; Correction of sst on OAFLUX grid 10 13 ; 11 ; :file:`${PROJECT_OD}/erai_sst_ 19890101_20091231_oafluxgrid.nc`14 ; :file:`${PROJECT_OD}/erai_sst_{yyyymmdd}_{yyyymmdd}_oafluxgrid.nc` 12 15 ; containing 13 16 ; ++ … … 17 20 ; Corrected sst on OAFLUX grid 18 21 ; is written in 19 ; :file:`${PROJECT_OD}/TropFlux_sst_ 19890101_20091231.nc`22 ; :file:`${PROJECT_OD}/TropFlux_sst_{yyyymmdd}_{yyyymmdd}.nc` 20 23 ; if this file not already exists. 21 24 ; 22 25 ; This output file 23 ; :file:`${PROJECT_OD}/TropFlux_sst_ 19890101_20091231.nc`26 ; :file:`${PROJECT_OD}/TropFlux_sst_{yyyymmdd}_{yyyymmdd}.nc` 24 27 ; will be used by 25 28 ; :func:`cronin_gustiness_ncdf` 26 29 ; and 27 ; : ref:`TropFlux_19890101_20091231.pro`.30 ; :func:`tropflux`. 28 31 ; 29 32 ; This output file 30 ; :file:`${PROJECT_OD}/TropFlux_sst_ 19890101_20091231.nc`33 ; :file:`${PROJECT_OD}/TropFlux_sst_{yyyymmdd}_{yyyymmdd}.nc` 31 34 ; is one of the deliverable of TropFlux. 32 35 ; … … 41 44 ; digraph sst_correction_ncdf { 42 45 ; 43 ; file_in [shape=ellipse,fontname=Courier,label="${PROJECT_OD}/erai_sst_ 19890101_20091231_oafluxgrid.nc"];44 ; file_out [shape=ellipse,fontname=Courier,label="${PROJECT_OD}/TropFlux_sst_ 19890101_20091231.nc"];46 ; file_in [shape=ellipse,fontname=Courier,label="${PROJECT_OD}/erai_sst_{yyyymmdd}_{yyyymmdd}_oafluxgrid.nc"]; 47 ; file_out [shape=ellipse,fontname=Courier,label="${PROJECT_OD}/TropFlux_sst_{yyyymmdd}_{yyyymmdd}.nc"]; 45 48 ; 46 49 ; sst_correction_ncdf [shape=box, … … 57 60 ; ======== 58 61 ; 62 ; :ref:`mooring_corrections` 63 ; 64 ; Used by :ref:`tropflux.sh` 65 ; 59 66 ; :ref:`project_profile.sh` 60 ;61 ; :ref:`mooring_corrections`62 67 ; 63 68 ; :func:`interp_erai_sst` … … 74 79 ; :: 75 80 ; 76 ; IDL> sst_correction_ncdf 81 ; IDL> yyyymmddb = 20000101L 82 ; IDL> yyyymmdde = 20001231L 83 ; IDL> result = sst_correction_ncdf(yyyymmddb, yyyymmdde) 84 ; IDL> print, result 77 85 ; 78 86 ; TODO … … 80 88 ; 81 89 ; make it work on loholt1 idl8: 82 ; reach end wh tiout error but with bad values (9.96921e+36 on time and sst, nb timestep 571690 ; reach end whithout error but with bad values (9.96921e+36 on time and sst, nb timestep 5716 83 91 ; 84 92 ; see commands:: … … 156 164 ; $URL$ 157 165 ; 166 ; - fplod 20120322 167 ; 168 ; * taking project_overwite into account 169 ; * try to add compile_opt seems to be incompatible with ncdf_quickwrite 170 ; * pro -> function 171 ; * hard coded da1 and da2 replaced by yyyymmddb and yyyymmdde parameters 172 ; * get rid of timegen 173 ; 158 174 ; - fplod 20110808T124236Z cratos (Linux) 159 175 ; … … 175 191 ; 176 192 ;- 177 pro sst_correction_ncdf 193 function sst_correction_ncdf $ 194 , yyyymmddb $ 195 , yyyymmdde 196 ; 197 ;++compile_opt idl2, strictarrsubs, logical_predicate 178 198 ; 179 199 @cm_4cal … … 183 203 @cm_project 184 204 ; 205 ; Return to caller if errors 206 ON_ERROR, 2 207 ; 208 result = -1 209 ; 210 usage = 'result = sst_correction_ncdf(yyyymmddb, yyyymmdde)' 211 nparam = N_PARAMS() 212 IF (nparam NE 2) THEN BEGIN 213 ras = report(['Incorrect number of arguments.' $ 214 + '!C' $ 215 + 'Usage : ' + usage]) 216 return, result 217 ENDIF 218 ; 185 219 ; test if ${PROJECT_OD} defined 186 220 CASE project_od_env OF … … 188 222 msg = 'eee : ${PROJECT_OD} is not defined' 189 223 ras = report(msg) 190 STOP224 return, result 191 225 END 192 226 ELSE: BEGIN … … 203 237 msg = 'eee : the directory' + iodirout + ' is not accessible.' 204 238 ras = report(msg) 205 STOP239 return, result 206 240 ENDIF 207 241 ; … … 210 244 msg = 'eee : the directory' + iodirout + ' was not found.' 211 245 ras = report(msg) 212 STOP 213 ENDIF 214 ; 215 da1=19880101 216 da2=20091231 246 return, result 247 ENDIF 217 248 ; 218 249 ; build data filename 219 filename='erai_sst_ 19890101_20091231_oafluxgrid.nc'250 filename='erai_sst_' + string(yyyymmddb,format='(I8.8)') + '_' + string(yyyymmdde,format='(I8.8)') + '_oafluxgrid.nc' 220 251 ; 221 252 ; check if this file exists … … 226 257 msg = 'eee : the file ' + fullfilename + ' was not found.' 227 258 ras = report(msg) 228 STOP259 return, result 229 260 ENDIF 230 261 ; 231 262 ; build output filename 232 filename_out = 'TropFlux_sst_ 19890101_20091231.nc'263 filename_out = 'TropFlux_sst_' + string(yyyymmddb,format='(I8.8)') + '_' + string(yyyymmdde,format='(I8.8)') + '.nc' 233 264 fullfilename_out = iodirout + filename_out 234 265 ; in order to avoid unexpected overwritten 235 IF ( FILE_TEST(fullfilename_out) EQ 1) THEN BEGIN266 IF ((FILE_TEST(fullfilename_out) EQ 1) AND (project_overwrite EQ 0)) THEN BEGIN 236 267 msg = 'eee : the file ' + fullfilename_out + ' already exists.' 237 268 ras = report(msg) 238 STOP269 return, result 239 270 ENDIF 240 271 ; 241 272 initncdf, fullfilename 242 sst=read_ncdf('sst',da1,da2,file=fullfilename,/nostr) 273 sst=read_ncdf('sst',yyyymmddb-.5d,yyyymmdde,file=fullfilename,/nostr) 274 timein=24.d*(time-julday(1,1,1957,0,0,0)) 275 jpt=n_elements(timein) 276 da=jul2date(time[0]) 277 cda0=string(da,format='(i8.8)') 278 da=jul2date(time[jpt-1]) 279 cda1=string(da,format='(i8.8)') 280 print, 'sst in sst_correction_ncdf first date ', cda0 281 print, 'sst in sst_correction_ncdf last date ' , cda1 282 ; 243 283 sst=sst-273.15 244 284 help, sst … … 250 290 caldat, time,mon,day,yea 251 291 ; 252 ; debug to understand some time value 9.9692100e+36 on idl 6 lohloht1253 print, time[7080:7090]292 ;++ debug to understand some time value 9.9692100e+36 on idl 6 lohloht1 293 ; print, time[7080:7090] when working on [19890101,19891231] 254 294 ; 255 295 sst_m=sst*0. 256 296 ; 257 297 for jt=0,jpt-1 do begin 258 jtt=(time (jt)-julday(1,1,yea(jt))) < 364259 t=reform(sst_mean (*,*))260 sst_m (*,*,jt)=t298 jtt=(time[jt]-julday(1,1,yea[jt])) < 364 299 t=reform(sst_mean[*,*]) 300 sst_m[*,*,jt]=t 261 301 endfor 262 302 help, sst_m … … 278 318 ; 279 319 ;writing field 280 lat=reform(gphit(0,0:jpj-1)) 281 lon=reform(glamt(0:jpi-1,0)) 282 time=timegen(7670, start=julday(1,1,1989,0), units='days') 283 jpt=n_elements(time) 284 ; 285 cda0=string(jul2date(time(0)),format='(i8.8)') 286 cda1=string(jul2date(time(jpt-1)),format='(i8.8)') 287 ; 288 time=time-julday(1,1,1950) 289 jpt=n_elements(time) 320 lat=reform(gphit[0,0:jpj-1]) 321 lon=reform(glamt[0:jpi-1,0]) 290 322 ; 291 323 ncfile='!' + fullfilename_out … … 296 328 sst_attr={units:'degK',missing_value:1.e20,long_name:'Sea Surface Temperature',short_name:'sst',axis:'TYX'} 297 329 ; 298 ncfields = 'sst[longitude,latitude, time]=sst_new:sst_attr; ' $330 ncfields = 'sst[longitude,latitude,*time]=sst_new:sst_attr; ' $ 299 331 + 'longitude[]=lon:lon_attr; ' $ 300 332 + 'latitude[]=lat:lat_attr; ' $ 301 + 'time[ *time]=time:time_attr ' $333 + 'time[]=timein:time_attr ' $ 302 334 + ' @ globattr' 303 335 ; 304 336 @ncdf_quickwrite 305 337 ; 338 result = 0 339 return, result 340 ; 306 341 end
Note: See TracChangeset
for help on using the changeset viewer.