From 922fb65f00af1fec0a2b54552e105f2f4bcec14c Mon Sep 17 00:00:00 2001
From: awalther <andi.walther@ssec.wisc.edu>
Date: Sun, 28 Oct 2018 09:06:38 -0500
Subject: [PATCH] add files
---
cws2__add_c.pro | 243 +++++++++++++++++++++++++++++++++++++++++++++
cws2__create_c.pro | 37 +++++++
cws2__define.pro | 243 ---------------------------------------------
3 files changed, 280 insertions(+), 243 deletions(-)
create mode 100644 cws2__add_c.pro
create mode 100644 cws2__create_c.pro
diff --git a/cws2__add_c.pro b/cws2__add_c.pro
new file mode 100644
index 0000000..8135c68
--- /dev/null
+++ b/cws2__add_c.pro
@@ -0,0 +1,243 @@
+;+
+;
+;-
+pro cws2::add_c, overpass = overpass, cth_reference = cth_reference, affi = affi, $
+ sensor = sensor, product=product, day = day, month = month, $
+ lon_SEVIRI=lon_SEVIRI, lat_SEVIRI=lat_SEVIRI
+
+ print, '*** START add_c *******!!!! '
+
+ if n_elements(affi) eq 0 then affi = self.affi
+ if n_elements(sensor) eq 0 then sensor = self.sensor
+ if n_elements(product) eq 0 then product = self.product
+
+ ; convert orbit number to starting time of the Atrain track
+ if n_elements(day) eq 0 or n_elements(month) eq 0 then begin
+ orbit_2_time, overpass, year, month, day, hour, minute
+ print, '... add data to AVAC-S level C file for '+affi+'/'+sensor+'/'+product+' ', year, month, day
+ endif
+ ;default,day,13
+ ;default,month,6
+
+ ; read time from level B data
+ bTime = self->get_data(product='YDT',/data)
+ ; time is given in seconds since 05th May 1968
+ time_ref = JULDAY(05,24,1968,0,0,0)
+ ;sec per day ; reference time
+ bTimeStart = bTime[0] / (3600D*24D) + time_ref
+ CALDAT, bTimeStart, Month, Day, Year, Hour, Minute, Second
+ print, '*** AVAC-S level B start date: ', Year, Month, Day, Hour, Minute, Second
+ bTimeEnd = bTime[-1] / (3600D*24D) + time_ref
+ CALDAT, bTimeEnd, Month, Day, Year, Hour, Minute, Second
+ print, '*** AVAC-S level B end date: ', Year, Month, Day, Hour, Minute, Second
+
+; ; seconds past midnight (of the particular day)
+; secDay = bTime mod (3600D*24D)
+;
+; srcSce_0 = floor(min(secDay[where( secDay ne !atrain.missing)])/60./15.)
+; srcSce_1 = ceil (max(secDay[where( secDay ne !atrain.missing)])/60./15.)
+;
+; print, min(secDay) / 3600.
+; print, max(secDay) / 3600.
+; print, ' *******!!!! ******* source scene:', srcSce_0, srcSce_1
+
+ ; start date
+ ; calculate seconds after a whole (quater of an hour)
+ ; = 15*60 seconds = 900 seconds
+ first_SEVIRI_date = bTime[0] - (bTime[0] mod (900D))
+ fSd = first_SEVIRI_date / (3600D*24D) + time_ref
+ CALDAT, fSd, Month, Day, Year, Hour, Minute, Second
+ print, '*** first SEVIRI date: ', Year, Month, Day, Hour, Minute, Second
+ last_SEVIRI_date = bTime[-1] + (900D - (bTime[-1] mod (900D)))
+ lSd = last_SEVIRI_date / (3600D*24D) + time_ref
+ CALDAT, lSd, Month, Day, Year, Hour, Minute, Second
+ print, '*** last SEVIRI date: ', Year, Month, Day, Hour, Minute, Second
+ n_disks = round ((last_SEVIRI_date - first_SEVIRI_date) / 900D)+1
+ print, '*** number of SEVIRI disks ', n_disks
+
+ ; delete variables
+ delvarx, cXele
+ delvarx, cXlin
+
+ ; read position of the ATRAIN PATH as index numbers of the MSG grid
+ ; from level B file
+ ; xele Column index of grid cell on SEVIRI disk (x)
+ ; xlin Row index of grid cell on SEVIRI disk (y)
+ cXele = self->_get_data(product='XELE',/gen)
+ cXlin = self->_get_data(product='XLIN',/gen)
+
+; help, cXele
+; help, cXlin
+; stop
+
+; !!! copied these line from UC_VDGEN/atrain_data__add_seviri.pro
+; might be useful for taking care of the parallax coorection
+;
+; ; ++++++++++++++++++++++++++++++++++++++++
+; ; get parallax indices for selected value of cth_Reference
+; ; ++++++++++++++++++++++++++++++++++++++++
+; self->get_parallax_indices, XLIN, XELE, XLIN_AREA, XELE_AREA, bCTH, cth_reference=cth_reference
+
+ ; define arrays: timeDiff compares time shift from different files
+; cTimeDiff = make_array(size(bTime,/dimension),/float,Value=999999.)
+ cTimeDiff = 450D ; 450s is equal to 7.5 min
+ cTime = make_array(size(bTime,/dimension),/double,Value=!ATRAIN.MISSING)
+ ; create data array along atrain path and write no data values everywhere
+ cData = make_array(size(bTime,/dimension),/float,Value=-1.)
+
+
+ ; -------- read SEVIRI data ----------
+ ; ------------------------------------
+
+ ;; define object cws_read
+ oCWS1 = obj_new('cws_read')
+ oCWS1 -> set_product, strlowcase(product)
+ oCWS1 -> set_group, affi
+
+ ;; read longitude and latitude
+ lon_SEVIRI = oCWS1->lon()
+ lat_SEVIRI = oCWS1->lat()
+ ; limit lon and lat to atrain track
+ lon = lon_SEVIRI[cxele,cxlin]
+ lat = lat_SEVIRI[cxele,cxlin]
+
+ ; search name of level C file
+ cFile = self->get_cFile()
+
+ ; counter for data points
+ n_data = 0
+
+ for iFile = 0, n_disks-1 do begin
+; for iFile = 0, 0 do begin
+ ; SEVIRI time of observation in days
+ itime = fSd + iFile * (900D / (3600D*24D))
+ CALDAT, itime, Month, Day, Year, Hour, Minute, Second
+ print, '*** read SEVIRI file: ', Year, Month, Day, Hour, Minute
+
+ ; precise time of SEVIRI observation in seconds
+ ; start time of the SEVIRI observation
+ t0 = first_SEVIRI_date + iFile * 900D
+ ; take care of scanning time
+ ; Scan of the disk from South to North takes 12min 40sec = 760sec
+ sTime = t0 + cXlin * (760./3712.)
+
+ index = where(abs(sTime - bTime) lt cTimeDiff, ntime)
+ print, 'ntime: ', ntime
+
+ ;; verbose output just for time checking
+ ;for tt=0, ntime-1, 10 do begin
+ ; CALDAT, sTime[index[tt]] /(3600D*24D)+time_ref, Month, Day, Year, Hour, Minute, Second
+ ; print, '*** SEVIRI obs: ', Year, Month, Day, Hour, Minute, Second
+ ; CALDAT, bTime[index[tt]] /(3600D*24D)+time_ref, Month, Day, Year, Hour, Minute, Second
+ ; print, ' ATRAIN obs: ', Year, Month, Day, Hour, Minute, Second
+ ; print, ' delta time: ', (sTime[index[tt]] - bTime[index[tt]]) / 60., ' min'
+ ;endfor
+
+ if nTime eq 0 then continue
+
+ print, '... read data from SEVIRI file for atrain index ', index[0], ' to ', index[-1]
+
+ ; add data points
+ n_data = n_data + ntime
+
+ ; read SEVIRI data
+ oCWS1 -> set_date, Year, Month, Day, Hour, Minute
+ data = oCWS1->get_data()
+
+ ; limit data to atrain path
+ data = data[cxele,cxlin]
+
+ print, '--- MIN/MAX(data)', min(data), max(data), min(data[index]), max(data[index])
+ cData[index] = data[index]
+
+ endfor
+
+;; old version looped over wrong files
+;-------------------------------------
+; for iFile = srcSce_0, srcSce_1 do begin
+;
+; ; calculate time of SEVIRI observation
+; ; Scan of the disk from South to North lasts 12 min 40sec = 760sec
+; sTimeNominal = (60.*60.*(ifile/4))+(60.*15.*(ifile mod 4))
+; sTime = sTimeNominal + cXlin * (760./3712.)
+; day_offset = bTime[0] - (bTime[0] mod (3600D * 24D))
+; sTime += day_offset
+;
+; index = where(abs(sTime - bTime) lt cTimeDiff, ntime)
+;
+; print, '((((( nTime: ', nTime
+; if nTime eq 0 then continue
+;
+; ; set new time shift and time
+; cTimeDiff[index] = (ABS(sTime-bTime))[index]
+; cTime[index] = sTime[index]
+;
+; ; define object cws_read
+; oCWS1 = obj_new('cws_read')
+; oCWS1 -> set_product, strlowcase(product)
+; oCWS1 -> set_group, affi
+; oCWS1 -> set_date,2008,month,day,ifile/4,15*(ifile mod 4)
+; data = oCWS1->get_data()
+;
+; if n_elements(data) eq 1 then continue
+;
+; ; limit data to atrain path
+; data = data[cxele,cxlin]
+;
+; ; read longitude and latitude
+; lon = oCWS1->lon()
+; lat = oCWS1->lat()
+;
+; ; destroy cloud workshop 1 object (was just there to read data)
+; obj_destroy, oCWS1
+;
+; ; limit lon and lat to atrain track
+; lon = lon[cxele,cxlin]
+; lat = lat[cxele,cxlin]
+;
+; print, '**** INDEX: ', index[0], index[-1]
+;
+; ;; data = cws2_read_data(affi,self.cws_sensor,product,xlin=cXlin, xele = cXele, lon=lon,lat=lat)
+; cData[index] = data[index]
+;
+; endfor
+
+
+ ; destroy cloud workshop 1 object (was just there to read data)
+ obj_destroy, oCWS1
+
+ print, '*** N_DATA', n_data
+
+ ; check for no data
+ if n_data eq 0 then begin
+ print, '*** no data found for '+affi+'/'+sensor+'/'+product+' ', year, month, day
+ return
+ endif
+
+ ; check for all data is not a value
+ if min(cdata) eq max(cdata) then begin
+ print, '*** Warning min = max = ', min(cdata),' for '+$
+ affi+'/'+sensor+'/'+product+' ', year, month, day
+ return
+ endif
+
+ ; lon/lat is written several times to the file, but this check is not very nice
+ ;if product eq 'cmb' then begin
+ ; procedure h5add in /usr/people/hamann/TOOLS/AVACS/stable-2.1/UC_VDGEN/atrain_data__h5add.pro
+ print, 'LON, add_c in cws2__define.pro'
+ self->h5add,cfile,lon, group='/'+affi+'/'+self.cws_sensor+'/'+'LON'
+ print, 'LAT, add_c in cws2__define.pro'
+ self->h5add,cfile,lat, group='/'+affi+'/'+self.cws_sensor+'/'+'LAT'
+ ;endif
+
+
+ ; add actual data to the level C file
+ print, '... add data to AVAC-S level C file ', ntime, ' datapoits '
+ print, '/'+affi+'/'+self.cws_sensor+'/'+product, ', add_c in cws2__define.pro'
+ self->h5add, cfile, cData, group='/'+affi+'/'+self.cws_sensor+'/'+product
+
+ print, ''
+ print, ''
+ print, ''
+
+end
diff --git a/cws2__create_c.pro b/cws2__create_c.pro
new file mode 100644
index 0000000..4f6cd66
--- /dev/null
+++ b/cws2__create_c.pro
@@ -0,0 +1,37 @@
+;+
+;
+;
+;
+;-
+pro cws2::create_c
+
+ self->set_sensor,'_GENERAL'
+ ; read number of x pixel in Seviri coordinates
+ xlin = self->get_data(product='XLIN',/data)
+ ; read number of y pixel in Seviri coordinates
+ xele = self->get_data(product='XELE',/data)
+ ; read time
+ time = self->get_data(product='YDT',/data)
+ ; get C-file
+ cFile = self->get_cfile()
+
+ self->h5add,cfile,xele $
+ ,Par_Name = 'Column on SEVIRI disk (x) ' $
+ ,sensor='_GENERAL' $
+ ,prd_code='XELE' $
+ ,Par_Info='Column index of grid cell on SEVIRI disk' $
+ ,/no_append
+
+ self->h5add,cfile,xlin $
+ ,Par_Name = 'Row on SEVIRI disk (y) ' $
+ ,sensor='_GENERAL' $
+ ,prd_code = 'XLIN' $
+ ,par_Info='Row index of grid cell on SEVIRI disk'
+
+ self->h5add,cfile,time $
+ ,sensor='_GENERAL' $
+ ,prd_code='YDT' $
+ ,par_Name = 'Year, day and time' $
+ ,par_Info='YDT given in seconds since midnight UT at the beginning of May 24 1968'
+
+end
diff --git a/cws2__define.pro b/cws2__define.pro
index c68e3c8..6aadccf 100644
--- a/cws2__define.pro
+++ b/cws2__define.pro
@@ -187,249 +187,6 @@ end
-;+
-;
-;-
-pro cws2::add_c, overpass = overpass, cth_reference = cth_reference, affi = affi, $
- sensor = sensor, product=product, day = day, month = month, $
- lon_SEVIRI=lon_SEVIRI, lat_SEVIRI=lat_SEVIRI
-
- print, '*** START add_c *******!!!! '
-
- if n_elements(affi) eq 0 then affi = self.affi
- if n_elements(sensor) eq 0 then sensor = self.sensor
- if n_elements(product) eq 0 then product = self.product
-
- ; convert orbit number to starting time of the Atrain track
- if n_elements(day) eq 0 or n_elements(month) eq 0 then begin
- orbit_2_time, overpass, year, month, day, hour, minute
- print, '... add data to AVAC-S level C file for '+affi+'/'+sensor+'/'+product+' ', year, month, day
- endif
- ;default,day,13
- ;default,month,6
-
- ; read time from level B data
- bTime = self->get_data(product='YDT',/data)
- ; time is given in seconds since 05th May 1968
- time_ref = JULDAY(05,24,1968,0,0,0)
- ;sec per day ; reference time
- bTimeStart = bTime[0] / (3600D*24D) + time_ref
- CALDAT, bTimeStart, Month, Day, Year, Hour, Minute, Second
- print, '*** AVAC-S level B start date: ', Year, Month, Day, Hour, Minute, Second
- bTimeEnd = bTime[-1] / (3600D*24D) + time_ref
- CALDAT, bTimeEnd, Month, Day, Year, Hour, Minute, Second
- print, '*** AVAC-S level B end date: ', Year, Month, Day, Hour, Minute, Second
-
-; ; seconds past midnight (of the particular day)
-; secDay = bTime mod (3600D*24D)
-;
-; srcSce_0 = floor(min(secDay[where( secDay ne !atrain.missing)])/60./15.)
-; srcSce_1 = ceil (max(secDay[where( secDay ne !atrain.missing)])/60./15.)
-;
-; print, min(secDay) / 3600.
-; print, max(secDay) / 3600.
-; print, ' *******!!!! ******* source scene:', srcSce_0, srcSce_1
-
- ; start date
- ; calculate seconds after a whole (quater of an hour)
- ; = 15*60 seconds = 900 seconds
- first_SEVIRI_date = bTime[0] - (bTime[0] mod (900D))
- fSd = first_SEVIRI_date / (3600D*24D) + time_ref
- CALDAT, fSd, Month, Day, Year, Hour, Minute, Second
- print, '*** first SEVIRI date: ', Year, Month, Day, Hour, Minute, Second
- last_SEVIRI_date = bTime[-1] + (900D - (bTime[-1] mod (900D)))
- lSd = last_SEVIRI_date / (3600D*24D) + time_ref
- CALDAT, lSd, Month, Day, Year, Hour, Minute, Second
- print, '*** last SEVIRI date: ', Year, Month, Day, Hour, Minute, Second
- n_disks = round ((last_SEVIRI_date - first_SEVIRI_date) / 900D)+1
- print, '*** number of SEVIRI disks ', n_disks
-
- ; delete variables
- delvarx, cXele
- delvarx, cXlin
-
- ; read position of the ATRAIN PATH as index numbers of the MSG grid
- ; from level B file
- ; xele Column index of grid cell on SEVIRI disk (x)
- ; xlin Row index of grid cell on SEVIRI disk (y)
- cXele = self->_get_data(product='XELE',/gen)
- cXlin = self->_get_data(product='XLIN',/gen)
-
-; help, cXele
-; help, cXlin
-; stop
-
-; !!! copied these line from UC_VDGEN/atrain_data__add_seviri.pro
-; might be useful for taking care of the parallax coorection
-;
-; ; ++++++++++++++++++++++++++++++++++++++++
-; ; get parallax indices for selected value of cth_Reference
-; ; ++++++++++++++++++++++++++++++++++++++++
-; self->get_parallax_indices, XLIN, XELE, XLIN_AREA, XELE_AREA, bCTH, cth_reference=cth_reference
-
- ; define arrays: timeDiff compares time shift from different files
-; cTimeDiff = make_array(size(bTime,/dimension),/float,Value=999999.)
- cTimeDiff = 450D ; 450s is equal to 7.5 min
- cTime = make_array(size(bTime,/dimension),/double,Value=!ATRAIN.MISSING)
- ; create data array along atrain path and write no data values everywhere
- cData = make_array(size(bTime,/dimension),/float,Value=-1.)
-
-
- ; -------- read SEVIRI data ----------
- ; ------------------------------------
-
- ;; define object cws_read
- oCWS1 = obj_new('cws_read')
- oCWS1 -> set_product, strlowcase(product)
- oCWS1 -> set_group, affi
-
- ;; read longitude and latitude
- lon_SEVIRI = oCWS1->lon()
- lat_SEVIRI = oCWS1->lat()
- ; limit lon and lat to atrain track
- lon = lon_SEVIRI[cxele,cxlin]
- lat = lat_SEVIRI[cxele,cxlin]
-
- ; search name of level C file
- cFile = self->get_cFile()
-
- ; counter for data points
- n_data = 0
-
- for iFile = 0, n_disks-1 do begin
-; for iFile = 0, 0 do begin
- ; SEVIRI time of observation in days
- itime = fSd + iFile * (900D / (3600D*24D))
- CALDAT, itime, Month, Day, Year, Hour, Minute, Second
- print, '*** read SEVIRI file: ', Year, Month, Day, Hour, Minute
-
- ; precise time of SEVIRI observation in seconds
- ; start time of the SEVIRI observation
- t0 = first_SEVIRI_date + iFile * 900D
- ; take care of scanning time
- ; Scan of the disk from South to North takes 12min 40sec = 760sec
- sTime = t0 + cXlin * (760./3712.)
-
- index = where(abs(sTime - bTime) lt cTimeDiff, ntime)
- print, 'ntime: ', ntime
-
- ;; verbose output just for time checking
- ;for tt=0, ntime-1, 10 do begin
- ; CALDAT, sTime[index[tt]] /(3600D*24D)+time_ref, Month, Day, Year, Hour, Minute, Second
- ; print, '*** SEVIRI obs: ', Year, Month, Day, Hour, Minute, Second
- ; CALDAT, bTime[index[tt]] /(3600D*24D)+time_ref, Month, Day, Year, Hour, Minute, Second
- ; print, ' ATRAIN obs: ', Year, Month, Day, Hour, Minute, Second
- ; print, ' delta time: ', (sTime[index[tt]] - bTime[index[tt]]) / 60., ' min'
- ;endfor
-
- if nTime eq 0 then continue
-
- print, '... read data from SEVIRI file for atrain index ', index[0], ' to ', index[-1]
-
- ; add data points
- n_data = n_data + ntime
-
- ; read SEVIRI data
- oCWS1 -> set_date, Year, Month, Day, Hour, Minute
- data = oCWS1->get_data()
-
- ; limit data to atrain path
- data = data[cxele,cxlin]
-
- print, '--- MIN/MAX(data)', min(data), max(data), min(data[index]), max(data[index])
- cData[index] = data[index]
-
- endfor
-
-;; old version looped over wrong files
-;-------------------------------------
-; for iFile = srcSce_0, srcSce_1 do begin
-;
-; ; calculate time of SEVIRI observation
-; ; Scan of the disk from South to North lasts 12 min 40sec = 760sec
-; sTimeNominal = (60.*60.*(ifile/4))+(60.*15.*(ifile mod 4))
-; sTime = sTimeNominal + cXlin * (760./3712.)
-; day_offset = bTime[0] - (bTime[0] mod (3600D * 24D))
-; sTime += day_offset
-;
-; index = where(abs(sTime - bTime) lt cTimeDiff, ntime)
-;
-; print, '((((( nTime: ', nTime
-; if nTime eq 0 then continue
-;
-; ; set new time shift and time
-; cTimeDiff[index] = (ABS(sTime-bTime))[index]
-; cTime[index] = sTime[index]
-;
-; ; define object cws_read
-; oCWS1 = obj_new('cws_read')
-; oCWS1 -> set_product, strlowcase(product)
-; oCWS1 -> set_group, affi
-; oCWS1 -> set_date,2008,month,day,ifile/4,15*(ifile mod 4)
-; data = oCWS1->get_data()
-;
-; if n_elements(data) eq 1 then continue
-;
-; ; limit data to atrain path
-; data = data[cxele,cxlin]
-;
-; ; read longitude and latitude
-; lon = oCWS1->lon()
-; lat = oCWS1->lat()
-;
-; ; destroy cloud workshop 1 object (was just there to read data)
-; obj_destroy, oCWS1
-;
-; ; limit lon and lat to atrain track
-; lon = lon[cxele,cxlin]
-; lat = lat[cxele,cxlin]
-;
-; print, '**** INDEX: ', index[0], index[-1]
-;
-; ;; data = cws2_read_data(affi,self.cws_sensor,product,xlin=cXlin, xele = cXele, lon=lon,lat=lat)
-; cData[index] = data[index]
-;
-; endfor
-
-
- ; destroy cloud workshop 1 object (was just there to read data)
- obj_destroy, oCWS1
-
- print, '*** N_DATA', n_data
-
- ; check for no data
- if n_data eq 0 then begin
- print, '*** no data found for '+affi+'/'+sensor+'/'+product+' ', year, month, day
- return
- endif
-
- ; check for all data is not a value
- if min(cdata) eq max(cdata) then begin
- print, '*** Warning min = max = ', min(cdata),' for '+$
- affi+'/'+sensor+'/'+product+' ', year, month, day
- return
- endif
-
- ; lon/lat is written several times to the file, but this check is not very nice
- ;if product eq 'cmb' then begin
- ; procedure h5add in /usr/people/hamann/TOOLS/AVACS/stable-2.1/UC_VDGEN/atrain_data__h5add.pro
- print, 'LON, add_c in cws2__define.pro'
- self->h5add,cfile,lon, group='/'+affi+'/'+self.cws_sensor+'/'+'LON'
- print, 'LAT, add_c in cws2__define.pro'
- self->h5add,cfile,lat, group='/'+affi+'/'+self.cws_sensor+'/'+'LAT'
- ;endif
-
-
- ; add actual data to the level C file
- print, '... add data to AVAC-S level C file ', ntime, ' datapoits '
- print, '/'+affi+'/'+self.cws_sensor+'/'+product, ', add_c in cws2__define.pro'
- self->h5add, cfile, cData, group='/'+affi+'/'+self.cws_sensor+'/'+product
-
- print, ''
- print, ''
- print, ''
-
-end
;+
--
GitLab