diff --git a/cws_read__define.pro b/cws_read__define.pro index a71581c67257b56f9849764616de4398e4e354ad..5be95ec7cba2f26e07c53d83a145c6b274fdcdbd 100755 --- a/cws_read__define.pro +++ b/cws_read__define.pro @@ -481,153 +481,6 @@ function cws_read::get_data, group=group, product=product, ice=ice, surf=surf, p end -;-------------------------------------------------------------------------------------------------------------------------------------- -pro cws_read::save_product, img, name, img_type, img_quality, $ - unit, accuracy, slope, offset, nodata_value_raw, nodata_value_phy, $ - version, used_data, remarks, $ - scale = scale, x_0 = x_0, y_0 = y_0, x_1 = x_1, y_1 = y_1, $ - irregular_path = irregular_path, irregular_filename = irregular_filename - - COMMON FEEDBACK, quiet, verbose, debug - - if verbose then print, '... start cws_read::save_product (cws_read__define.pro)' - - ;image dimensionen pruefen - dum_size = size(img) - n_dims = dum_size[0] - case n_dims of - 2: begin - true = 0 - end - 3: begin - idx = where(dum_size[1:3] eq 3, anz) - if anz ne 1 then begin - print, '*** Error in save_product: wrong image dimensions (product: "' + name + '") -> skip' - self -> print_log, '*** Error in save_product: wrong image dimensions (product: "' + name + '") -> skip', log_mode = 32 - return - endif - true = 1 + idx[0] - end - else: begin - print, '*** Error in save_product: image dimension has to be 2 or 3 (product: "' + name + '") -> skip' - self -> print_log, '*** Error in save_product: image dimension has to be 2 or 3 (product: "' + name + '") -> skip', log_mode = 32 - return - endelse - endcase - - ; calculate/read limits for scale and img - scale = advanced_keyword_set(scale) ? scale : self.scale - x_0 = advanced_keyword_set(x_0) ? x_0 : self.x_0 - y_0 = advanced_keyword_set(y_0) ? y_0 : self.y_0 - x_1 = advanced_keyword_set(x_1) ? x_1 : self.x_1 - y_1 = advanced_keyword_set(y_1) ? y_1 : self.y_1 - dum_dim = dum_size(1:dum_size[0]) - dum_dim = dum_dim[where(dum_dim gt 3)] - if scale eq 0 then scale = 1.0 - if x_0 eq 0 and x_1 eq 0 then x_1 = dum_dim[0] - 1 - if y_0 eq 0 and y_1 eq 0 then y_1 = dum_dim[0] - 1 - ;pfad ueberpruefen und ggf. erstellen - path = keyword_set(irregular_path) ? irregular_path : (self -> build_product_path()) - - if not file_test(path,/DIR) then begin - file_mkdir, path - ;if not ok then begin & self -> print_log, 'path creation failed: ' + path + ' -> skip', log_mode = 32 & return & endif - endif - - ; create filenames - img_file = keyword_set(irregular_filename) ? (irregular_filename + '.' + img_type) : (self -> build_product_filename(name, img_type)) - xml_file = keyword_set(irregular_filename) ? (irregular_filename + '.xml') : (self -> build_product_filename(name, 'xml')) - - catch, error - if error eq 0 then begin - ; write img_file - if not quiet then print, ' write image file: ', img_file - ;window, 11 - ;view2d,img,/cool,/colo,no_data_idx=where(img le 0) - case img_type of - 'jpg': begin - self -> print_log, img_file, log_mode = 16 - write_jpeg, path + img_file, img, quality = img_quality, true = true - depth = 8 - end - 'png': begin - self -> print_log, img_file, log_mode = 16 - ;help, img - write_png, path + img_file, img - print, path + img_file - dum = size(img, /type) - depth = ((dum eq 2) or (dum eq 3) or ((dum ge 12) and (dum le 15))) ? 16 : 8 - ;help, img - ;view2d, img, /cool, /colo ;, no_data_idx=where(img le 0) - ;stop - end - else : begin - print, '*** Error in save_product: not yet supported image type: ' + type + ' -> skip' - self -> print_log, 'not yet supported image type: ' + type + ' -> skip', log_mode = 32 - return - end - endcase - - ;xml_file schreiben - if not quiet then print, ' write xml file: ', xml_file - if n_elements(used_data) eq 0 then print, '*** Error in save_product: no information about used_data given' - openw, lun_1, path + xml_file, /get_lun - printf, lun_1, "<?xml version='1.0' standalone='yes'?>" - print, '<MSGProduct Name="' + name - print, '" File="' + img_file - print, '" DateOfProcessing="' + systime() - print, '" ProcessorVersion="' + version - print, '" UsedData="' + used_data - print, '" Remarks="' + remarks + '">' - printf, lun_1, '<MSGProduct Name="' + name + $ - '" File="' + img_file + $ - '" DateOfProcessing="' + systime() + $ - '" ProcessorVersion="' + version + $ - '" UsedData="' + used_data + $ - '" Remarks="' + remarks + '">' - printf, lun_1, ' <ImageInformation Type="' + img_type + $ - '" BitPerPixel="' + strcompress(string(depth, format = '(i2)'), /remove_all) + $ - '" ImgChannels="' + ((true ne 0) ? '3' : '1') + $ - '" Quality="' + strcompress(string(img_quality, format = '(f15.5)'), /remove_all) + $ - '" Scale="' + strcompress(string(scale, format = '(f15.5)'), /remove_all) + $ - '" X0="' + strcompress(string(x_0, format = '(i5)'), /remove_all) + $ - '" Y0="' + strcompress(string(y_0, format = '(i5)'), /remove_all) + $ - '" X1="' + strcompress(string(x_1, format = '(i5)'), /remove_all) + $ - '" Y1="' + strcompress(string(y_1, format = '(i5)'), /remove_all) + $ - '">' - print, ' </ImageInformation>' - print, ' <DataRepresentation Slope="' + strcompress(string(slope, format = '(f15.5)'), /remove_all) + $ - '" Offset="' + strcompress(string(offset, format = '(f15.5)'), /remove_all) + $ - '" NoDataValueRaw="' + strcompress(string(nodata_value_raw, format = '(f15.5)'), /remove_all) + $ - '" NoDataValuePhysical="' + strcompress(string(nodata_value_phy, format = '(f15.5)'), /remove_all) + $ - '" PhysicalUnit="' + unit + $ - '" Accuracy="' + strcompress(string(accuracy, format = '(f15.5)'), /remove_all) + $ - '">' - print, ' </DataRepresentation>' - print, '</MSGProduct>' - printf, lun_1, ' </DataRepresentation>' - printf, lun_1, '</MSGProduct>' - printf, lun_1, ' </ImageInformation>' - printf, lun_1, ' <DataRepresentation Slope="' + strcompress(string(slope, format = '(f15.5)'), /remove_all) + $ - '" Offset="' + strcompress(string(offset, format = '(f15.5)'), /remove_all) + $ - '" NoDataValueRaw="' + strcompress(string(nodata_value_raw, format = '(f15.5)'), /remove_all) + $ - '" NoDataValuePhysical="' + strcompress(string(nodata_value_phy, format = '(f15.5)'), /remove_all) + $ - '" PhysicalUnit="' + unit + $ - '" Accuracy="' + strcompress(string(accuracy, format = '(f15.5)'), /remove_all) + $ - '">' - printf, lun_1, ' </DataRepresentation>' - printf, lun_1, '</MSGProduct>' - free_lun, lun_1 - - endif else begin - print, '*** Error in save_product: can not write product to harddrive: ' + path + img_file + ' -> skip' - self -> print_log, '*** Error in save_product:can not write product to harddrive: ' + path + img_file + ' -> skip', log_mode = 32 - return - endelse - catch, /cancel - -end -;-------------------------------------------------------------------------------------------------------------------------------------- ;-------------------------------------------------------------------------------------------------------------------------------------- diff --git a/cws_read__readData_TPS.pro b/cws_read__readData_TPS.pro new file mode 100644 index 0000000000000000000000000000000000000000..6578a51645f0ae9fef17dd68e73165ab0d903cfd --- /dev/null +++ b/cws_read__readData_TPS.pro @@ -0,0 +1,31 @@ +; reads TROPOS DATA +; +function cws_read::readData_TPS + CASE strLowCase(self.product) OF + 'cmb': begin + self.infile = FILE_SEARCH(self->build_source_path() $ + + '/TROPOS_MSI_CM_*'+ self->get_date(/string,/KMR) +'*.h5',c=c) + + fileId = h5f_open(self.infile) + dataId = h5d_open(fileID,'ScienceData/cloud_mask') + data = h5d_read(dataID) + img = make_array(3712,3712,/byte,value=0.) + img = (data eq 0 or data eq 1) + 100 * ( data eq 2 or data eq 3) + img = rotate(img,7) + + return,img + end + + 'ctp': begin + + + + + end + + + endcase + + + +end function diff --git a/cws_read__readdata.pro b/cws_read__readdata.pro index c81a6f8aa7077e80a31df956967cf5051f0a0b43..2141ad90a6bd5a0e654951d489fa8b0c1515cc71 100644 --- a/cws_read__readdata.pro +++ b/cws_read__readdata.pro @@ -28,50 +28,11 @@ FUNCTION cws_read::readData, version=version CASE self.grp OF 'TPS':begin - CASE strLowCase(self.product) OF - 'cmb': begin - self.infile = FILE_SEARCH(self->build_source_path() $ - + '/TROPOS_MSI_CM_*'+ self->get_date(/string,/KMR) +'*.h5',c=c) - - fileId = h5f_open(self.infile) - dataId = h5d_open(fileID,'ScienceData/cloud_mask') - data = h5d_read(dataID) - img = make_array(3712,3712,/byte,value=0.) - img = (data eq 0 or data eq 1) + 100 * ( data eq 2 or data eq 3) - img = rotate(img,7) - - return,img - end - - - endcase + img = self.readData_TPS() end 'KNM':BEGIN - CASE strLowCase(self.product) OF - 'cmb': begin - self.infile = FILE_SEARCH(self->build_source_path() $ - + '/msg2_aodc*'+ self->get_date(/string,/KMR) +'*064-163.nc',c=c) - if c ne 1 then begin - print,'no source file: ', self.build_source_path() +'/msg2_aodc*'+self.get_date(/string,/KMR) $ - + '*064-163.nc' - return,-1 - endif - NCDF_GET, self.infile,'cldmask',cld - cmb = cld['cldmask','value'] - img = make_array(3712,3712,/byte,value=0.) - img = (cmb eq 0) + 100 * ( cmb eq 1) - img = rotate(img,7) - - return,img - - - end - - 'cod': - - endcase - + img = self.readdata_KNM() END @@ -435,8 +396,8 @@ FUNCTION cws_read::readData, version=version data = read_binary(self.infile,data_dims=[3712,3712],data_type=2,endian='big') data = rotate(data,7) - cmb = self->get_data(pr='cmb') - self.product='cph' + cmb = self->get_data(pr='cmb') + self.product='cph' ;window,1 ;view2d, data, /cool, /colo, no_data_idx=where(data lt 0) @@ -874,189 +835,8 @@ FUNCTION cws_read::readData, version=version ; 'DLR':BEGIN - if keyword_set(uncertainty) then error_message_uncertainty, self.grp - - CASE self.product of - - 'cmb':BEGIN - - self.infile = file_search(self->build_source_path()+self->get_date(/string)+'_MSG2_HRIT_CLM___.dat',c=c) - IF c eq 0 THEN BEGIN - IF verbose THEN BEGIN - PRINT, '... Error in cws_read__readdata' - PRINT, ' can not read ', self->build_source_path()+self->get_date(/string)+'_MSG2_HRIT_CLM___.dat' - ENDIF - RETURN, -1 - ENDIF ELSE IF verbose THEN PRINT, '... read ', self.infile - - data = read_binary(self.infile,data_dims=[3712,3712],data_type=1) - data = rotate(temporary(data),5) - ncdf_get_field, '/usr/people/hamann/TOOLS/cws3/trunk/MSG_data/MSG_disc.nc', 'disc', disc - ; [0-noData - img = 1*(disc eq 1) $ ; 1-noCloud - + 2*(data gt 80) ; 2-cloud] - RETURN, img - END - - 'ctp': BEGIN - self.infile = file_search(self->build_source_path()+self->get_date(/string)+'_MSG2_HRIT_CTH___.dat',c=c) - IF c eq 0 THEN BEGIN - IF verbose THEN BEGIN - PRINT, '... Error in cws_read__readdata' - PRINT, ' can not read ', self->build_source_path()+self->get_date(/string)+'_MSG2_HRIT_CTH___.dat' - ENDIF - RETURN, -1 - ENDIF ELSE IF verbose THEN PRINT, '... read ', self.infile - - data = read_binary(self.infile,data_dims=[3712,3712],data_type=2) - data = rotate(temporary(data),5) - - idx = where(data lt 100,c) - IF c gt 0 THEN data[idx] = -1 - - RETURN,data - end - - 'cth': BEGIN - self.infile = file_search(self->build_source_path()+self->get_date(/string)+'_MSG2_HRIT_CTZ___.dat',c=c) - IF c eq 0 THEN BEGIN - IF verbose THEN BEGIN - PRINT, '... Error in cws_read__readdata' - PRINT, ' can not read ', self->build_source_path()+self->get_date(/string)+'_MSG2_HRIT_CTZ___.dat' - ENDIF - RETURN, -1 - ENDIF ELSE IF verbose THEN PRINT, '... read ', self.infile - - data = read_binary(self.infile,data_dims=[3712,3712],data_type=2) - data = rotate(temporary(data),5) - - idx = where(data lt 100,c) - IF c gt 0 THEN data[idx] = -1 - - RETURN,data - end - - 'ctt' : BEGIN - self.infile = file_search(self->build_source_path()+self->get_date(/string)+'_MSG2_HRIT_CTT___.dat',c=c) - IF c eq 0 THEN BEGIN - IF verbose THEN BEGIN - PRINT, '... Error in cws_read__readdata' - PRINT, ' can not read ', self->build_source_path()+self->get_date(/string)+'_MSG2_HRIT_CTT___.dat' - ENDIF - RETURN, -1 - ENDIF ELSE IF verbose THEN PRINT, '... read ', self.infile - - data = read_binary(self.infile,data_dims=[3712,3712],data_type=2) - data = rotate(temporary(data),5) - - idx = where(data eq 0,c) - img = temporary(data)/10. - IF c gt 0 THEN img[idx] = -1 - - RETURN,img - end - - 'cod': BEGIN - file_ic = file_search(self->build_source_path()+self->get_date(/string)+'_MSG2_HRIT_TAU___IC_.dat',c=c) - IF c eq 0 THEN BEGIN - IF verbose THEN BEGIN - PRINT, '... Error in cws_read__readdata' - PRINT, ' can not read ', self->build_source_path()+self->get_date(/string)+'_MSG2_HRIT_TAU___IC_.dat' - ENDIF - RETURN, -1 - ENDIF ELSE IF verbose THEN PRINT, '... read ', file_ic - data_ic = read_binary(file_ic,data_dims=[3712,3712],data_type=4) - - file_wc = file_search(self->build_source_path()+self->get_date(/string)+'_MSG2_HRIT_TAU___WC_.dat',c=c) - IF c eq 0 THEN BEGIN - IF verbose THEN BEGIN - PRINT, '... Error in cws_read__readdata' - PRINT, ' can not read ', self->build_source_path()+self->get_date(/string)+'_MSG2_HRIT_TAU___WC_.dat' - ENDIF - RETURN, -1 - ENDIF ELSE IF verbose THEN PRINT, '... read ', file_wc - data_wc = read_binary(file_wc,data_dims=[3712,3712],data_type=4) - - img = data_wc > data_ic - img = rotate(img,5) - - ; correct some non existing data - cmb = self->get_data(product='cmb',group='DLR') - self.product = 'cod' - nodata_idx = where(cmb lt 1.5, nodata_anz) - IF nodata_anz gt 0 THEN img[nodata_idx]=-1. - - RETURN, img - END - - 'ref': BEGIN - file_ic = file_search(self->build_source_path()+self->get_date(/string)+'_MSG2_HRIT_REFF__IC_.dat',c=c) - IF c eq 0 THEN BEGIN - IF verbose THEN BEGIN - PRINT, '... Error in cws_read__readdata' - PRINT, ' can not read ', self->build_source_path()+self->get_date(/string)+'_MSG2_HRIT_REFF__IC_.dat' - ENDIF - RETURN, -1 - ENDIF ELSE IF verbose THEN PRINT, '... read ', file_ic - data_ic = read_binary(file_ic,data_dims=[3712,3712],data_type=4) - - file_wc = file_search(self->build_source_path()+self->get_date(/string)+'_MSG2_HRIT_REFF__WC_.dat',c=c) - IF c eq 0 THEN BEGIN - IF verbose THEN BEGIN - PRINT, '... Error in cws_read__readdata' - PRINT, ' can not read ', self->build_source_path()+self->get_date(/string)+'_MSG2_HRIT_REFF__WC_.dat' - ENDIF - RETURN, -1 - ENDIF ELSE IF verbose THEN PRINT, '... read ', file_wc - data_wc = read_binary(file_wc,data_dims=[3712,3712],data_type=4) - - img = data_wc > data_ic - img = rotate(img,5) - - cmb = self->get_data(product='cmb',group='DLR') - self.product = 'ref' - - nodata_idx = where(cmb lt 1.5,nodata_anz) - IF nodata_anz gt 0 THEN img[nodata_idx]=-1. - RETURN,img - end - - 'cph': BEGIN - self.infile = file_search(self->build_source_path()+self->get_date(/string)+'_MSG2_HRIT_PHASE_.dat',c=c) - IF c eq 0 THEN BEGIN - IF verbose THEN BEGIN - PRINT, '... Error in cws_read__readdata' - PRINT, ' can not read ', self->build_source_path()+self->get_date(/string)+'_MSG2_HRIT_CTT___.dat' - ENDIF - RETURN, -1 - ENDIF ELSE IF verbose THEN PRINT, '... read ', self.infile - data = read_binary(self.infile,data_dims=[3712,3712],data_type=1) - data = rotate(data,5) - - ;print, 'DLR raw data: ', min(data), max(data) - ;window,1 - ;view2d, data, /cool, /colo ;, no_data_idx=where(data eq 0) - histo = histogram(data, LOCATIONS=loc) - for i = 0, n_elements(histo)-1 do print, loc[i], histo[i] - - img = -1.*(data eq 0) +0.*(data eq 1) +100.*(data eq 2) - ; outside disc ; water cloud ; ice cloud - - ;; why is that needed??? - ;; read cloud binary mask as help file - ;cmb = self->get_data(product='cmb',group='DLR') - ;cmb = round(cmb) - ;idx = where(cmb eq 1) - ;data[idx] = 4 - ;; set product back to cloud phase - ;self.product = 'cph' - - ;window, 3 - ;view2d, cmb, /cool, /color, /no_axes, /asp, /full_screen - ;window, 4 - ;view2d, data, /cool, /color, /no_axes, /asp, /full_screen - - RETURN,img + img = self.readdata_dlr() + END ELSE: BEGIN diff --git a/cws_read__readdata_DLR.pro b/cws_read__readdata_DLR.pro new file mode 100644 index 0000000000000000000000000000000000000000..fea57ffac567004a6e8f37b6178ea6ad424e7170 --- /dev/null +++ b/cws_read__readdata_DLR.pro @@ -0,0 +1,190 @@ +; reads DLR data version 2008 +; +function cws_read::readdata_DLR + +if keyword_set(uncertainty) then error_message_uncertainty, self.grp + + CASE self.product of + + 'cmb':BEGIN + + self.infile = file_search(self->build_source_path()+self->get_date(/string)+'_MSG2_HRIT_CLM___.dat',c=c) + IF c eq 0 THEN BEGIN + IF verbose THEN BEGIN + PRINT, '... Error in cws_read__readdata' + PRINT, ' can not read ', self->build_source_path()+self->get_date(/string)+'_MSG2_HRIT_CLM___.dat' + ENDIF + RETURN, -1 + ENDIF ELSE IF verbose THEN PRINT, '... read ', self.infile + + data = read_binary(self.infile,data_dims=[3712,3712],data_type=1) + data = rotate(temporary(data),5) + ncdf_get_field, '/usr/people/hamann/TOOLS/cws3/trunk/MSG_data/MSG_disc.nc', 'disc', disc + ; [0-noData + img = 1*(disc eq 1) $ ; 1-noCloud + + 2*(data gt 80) ; 2-cloud] + RETURN, img + END + + 'ctp': BEGIN + self.infile = file_search(self->build_source_path()+self->get_date(/string)+'_MSG2_HRIT_CTH___.dat',c=c) + IF c eq 0 THEN BEGIN + IF verbose THEN BEGIN + PRINT, '... Error in cws_read__readdata' + PRINT, ' can not read ', self->build_source_path()+self->get_date(/string)+'_MSG2_HRIT_CTH___.dat' + ENDIF + RETURN, -1 + ENDIF ELSE IF verbose THEN PRINT, '... read ', self.infile + + data = read_binary(self.infile,data_dims=[3712,3712],data_type=2) + data = rotate(temporary(data),5) + + idx = where(data lt 100,c) + IF c gt 0 THEN data[idx] = -1 + + RETURN,data + end + + 'cth': BEGIN + self.infile = file_search(self->build_source_path()+self->get_date(/string)+'_MSG2_HRIT_CTZ___.dat',c=c) + IF c eq 0 THEN BEGIN + IF verbose THEN BEGIN + PRINT, '... Error in cws_read__readdata' + PRINT, ' can not read ', self->build_source_path()+self->get_date(/string)+'_MSG2_HRIT_CTZ___.dat' + ENDIF + RETURN, -1 + ENDIF ELSE IF verbose THEN PRINT, '... read ', self.infile + + data = read_binary(self.infile,data_dims=[3712,3712],data_type=2) + data = rotate(temporary(data),5) + + idx = where(data lt 100,c) + IF c gt 0 THEN data[idx] = -1 + + RETURN,data + end + + 'ctt' : BEGIN + self.infile = file_search(self->build_source_path()+self->get_date(/string)+'_MSG2_HRIT_CTT___.dat',c=c) + IF c eq 0 THEN BEGIN + IF verbose THEN BEGIN + PRINT, '... Error in cws_read__readdata' + PRINT, ' can not read ', self->build_source_path()+self->get_date(/string)+'_MSG2_HRIT_CTT___.dat' + ENDIF + RETURN, -1 + ENDIF ELSE IF verbose THEN PRINT, '... read ', self.infile + + data = read_binary(self.infile,data_dims=[3712,3712],data_type=2) + data = rotate(temporary(data),5) + + idx = where(data eq 0,c) + img = temporary(data)/10. + IF c gt 0 THEN img[idx] = -1 + + RETURN,img + end + + 'cod': BEGIN + file_ic = file_search(self->build_source_path()+self->get_date(/string)+'_MSG2_HRIT_TAU___IC_.dat',c=c) + IF c eq 0 THEN BEGIN + IF verbose THEN BEGIN + PRINT, '... Error in cws_read__readdata' + PRINT, ' can not read ', self->build_source_path()+self->get_date(/string)+'_MSG2_HRIT_TAU___IC_.dat' + ENDIF + RETURN, -1 + ENDIF ELSE IF verbose THEN PRINT, '... read ', file_ic + data_ic = read_binary(file_ic,data_dims=[3712,3712],data_type=4) + + file_wc = file_search(self->build_source_path()+self->get_date(/string)+'_MSG2_HRIT_TAU___WC_.dat',c=c) + IF c eq 0 THEN BEGIN + IF verbose THEN BEGIN + PRINT, '... Error in cws_read__readdata' + PRINT, ' can not read ', self->build_source_path()+self->get_date(/string)+'_MSG2_HRIT_TAU___WC_.dat' + ENDIF + RETURN, -1 + ENDIF ELSE IF verbose THEN PRINT, '... read ', file_wc + data_wc = read_binary(file_wc,data_dims=[3712,3712],data_type=4) + + img = data_wc > data_ic + img = rotate(img,5) + + ; correct some non existing data + cmb = self->get_data(product='cmb',group='DLR') + self.product = 'cod' + nodata_idx = where(cmb lt 1.5, nodata_anz) + IF nodata_anz gt 0 THEN img[nodata_idx]=-1. + + RETURN, img + END + + 'ref': BEGIN + file_ic = file_search(self->build_source_path()+self->get_date(/string)+'_MSG2_HRIT_REFF__IC_.dat',c=c) + IF c eq 0 THEN BEGIN + IF verbose THEN BEGIN + PRINT, '... Error in cws_read__readdata' + PRINT, ' can not read ', self->build_source_path()+self->get_date(/string)+'_MSG2_HRIT_REFF__IC_.dat' + ENDIF + RETURN, -1 + ENDIF ELSE IF verbose THEN PRINT, '... read ', file_ic + data_ic = read_binary(file_ic,data_dims=[3712,3712],data_type=4) + + file_wc = file_search(self->build_source_path()+self->get_date(/string)+'_MSG2_HRIT_REFF__WC_.dat',c=c) + IF c eq 0 THEN BEGIN + IF verbose THEN BEGIN + PRINT, '... Error in cws_read__readdata' + PRINT, ' can not read ', self->build_source_path()+self->get_date(/string)+'_MSG2_HRIT_REFF__WC_.dat' + ENDIF + RETURN, -1 + ENDIF ELSE IF verbose THEN PRINT, '... read ', file_wc + data_wc = read_binary(file_wc,data_dims=[3712,3712],data_type=4) + + img = data_wc > data_ic + img = rotate(img,5) + + cmb = self->get_data(product='cmb',group='DLR') + self.product = 'ref' + + nodata_idx = where(cmb lt 1.5,nodata_anz) + IF nodata_anz gt 0 THEN img[nodata_idx]=-1. + RETURN,img + end + + 'cph': BEGIN + self.infile = file_search(self->build_source_path()+self->get_date(/string)+'_MSG2_HRIT_PHASE_.dat',c=c) + IF c eq 0 THEN BEGIN + IF verbose THEN BEGIN + PRINT, '... Error in cws_read__readdata' + PRINT, ' can not read ', self->build_source_path()+self->get_date(/string)+'_MSG2_HRIT_CTT___.dat' + ENDIF + RETURN, -1 + ENDIF ELSE IF verbose THEN PRINT, '... read ', self.infile + data = read_binary(self.infile,data_dims=[3712,3712],data_type=1) + data = rotate(data,5) + + ;print, 'DLR raw data: ', min(data), max(data) + ;window,1 + ;view2d, data, /cool, /colo ;, no_data_idx=where(data eq 0) + histo = histogram(data, LOCATIONS=loc) + for i = 0, n_elements(histo)-1 do print, loc[i], histo[i] + + img = -1.*(data eq 0) +0.*(data eq 1) +100.*(data eq 2) + ; outside disc ; water cloud ; ice cloud + + ;; why is that needed??? + ;; read cloud binary mask as help file + ;cmb = self->get_data(product='cmb',group='DLR') + ;cmb = round(cmb) + ;idx = where(cmb eq 1) + ;data[idx] = 4 + ;; set product back to cloud phase + ;self.product = 'cph' + + ;window, 3 + ;view2d, cmb, /cool, /color, /no_axes, /asp, /full_screen + ;window, 4 + ;view2d, data, /cool, /color, /no_axes, /asp, /full_screen + + RETURN,img + + +end diff --git a/cws_read__readdata_KNM.pro b/cws_read__readdata_KNM.pro new file mode 100644 index 0000000000000000000000000000000000000000..cceef298fa3b19ff6ee0ef5896b250e860c2b9d6 --- /dev/null +++ b/cws_read__readdata_KNM.pro @@ -0,0 +1,27 @@ +; rreads ICWG2 KNMI data +; +function cws_read::readData_KNM +CASE strLowCase(self.product) OF + 'cmb': begin + self.infile = FILE_SEARCH(self->build_source_path() $ + + '/msg2_aodc*'+ self->get_date(/string,/KMR) +'*064-163.nc',c=c) + if c ne 1 then begin + print,'no source file: ', self.build_source_path() +'/msg2_aodc*'+self.get_date(/string,/KMR) $ + + '*064-163.nc' + return,-1 + endif + NCDF_GET, self.infile,'cldmask',cld + cmb = cld['cldmask','value'] + img = make_array(3712,3712,/byte,value=0.) + img = (cmb eq 0) + 100 * ( cmb eq 1) + img = rotate(img,7) + + return,img + + + end + + 'cod': + + endcase + end function diff --git a/cws_read__save_product.pro b/cws_read__save_product.pro new file mode 100644 index 0000000000000000000000000000000000000000..ff467b90519cd990916618e09d2b329920d0c9e8 --- /dev/null +++ b/cws_read__save_product.pro @@ -0,0 +1,147 @@ +;-------------------------------------------------------------------------------------------------------------------------------------- +pro cws_read::save_product, img, name, img_type, img_quality, $ + unit, accuracy, slope, offset, nodata_value_raw, nodata_value_phy, $ + version, used_data, remarks, $ + scale = scale, x_0 = x_0, y_0 = y_0, x_1 = x_1, y_1 = y_1, $ + irregular_path = irregular_path, irregular_filename = irregular_filename + + COMMON FEEDBACK, quiet, verbose, debug + + if verbose then print, '... start cws_read::save_product (cws_read__define.pro)' + + ;image dimensionen pruefen + dum_size = size(img) + n_dims = dum_size[0] + case n_dims of + 2: begin + true = 0 + end + 3: begin + idx = where(dum_size[1:3] eq 3, anz) + if anz ne 1 then begin + print, '*** Error in save_product: wrong image dimensions (product: "' + name + '") -> skip' + self -> print_log, '*** Error in save_product: wrong image dimensions (product: "' + name + '") -> skip', log_mode = 32 + return + endif + true = 1 + idx[0] + end + else: begin + print, '*** Error in save_product: image dimension has to be 2 or 3 (product: "' + name + '") -> skip' + self -> print_log, '*** Error in save_product: image dimension has to be 2 or 3 (product: "' + name + '") -> skip', log_mode = 32 + return + endelse + endcase + + ; calculate/read limits for scale and img + scale = advanced_keyword_set(scale) ? scale : self.scale + x_0 = advanced_keyword_set(x_0) ? x_0 : self.x_0 + y_0 = advanced_keyword_set(y_0) ? y_0 : self.y_0 + x_1 = advanced_keyword_set(x_1) ? x_1 : self.x_1 + y_1 = advanced_keyword_set(y_1) ? y_1 : self.y_1 + dum_dim = dum_size(1:dum_size[0]) + dum_dim = dum_dim[where(dum_dim gt 3)] + if scale eq 0 then scale = 1.0 + if x_0 eq 0 and x_1 eq 0 then x_1 = dum_dim[0] - 1 + if y_0 eq 0 and y_1 eq 0 then y_1 = dum_dim[0] - 1 + ;pfad ueberpruefen und ggf. erstellen + path = keyword_set(irregular_path) ? irregular_path : (self -> build_product_path()) + + if not file_test(path,/DIR) then begin + file_mkdir, path + ;if not ok then begin & self -> print_log, 'path creation failed: ' + path + ' -> skip', log_mode = 32 & return & endif + endif + + ; create filenames + img_file = keyword_set(irregular_filename) ? (irregular_filename + '.' + img_type) : (self -> build_product_filename(name, img_type)) + xml_file = keyword_set(irregular_filename) ? (irregular_filename + '.xml') : (self -> build_product_filename(name, 'xml')) + + catch, error + if error eq 0 then begin + ; write img_file + if not quiet then print, ' write image file: ', img_file + ;window, 11 + ;view2d,img,/cool,/colo,no_data_idx=where(img le 0) + case img_type of + 'jpg': begin + self -> print_log, img_file, log_mode = 16 + write_jpeg, path + img_file, img, quality = img_quality, true = true + depth = 8 + end + 'png': begin + self -> print_log, img_file, log_mode = 16 + ;help, img + write_png, path + img_file, img + print, path + img_file + dum = size(img, /type) + depth = ((dum eq 2) or (dum eq 3) or ((dum ge 12) and (dum le 15))) ? 16 : 8 + ;help, img + ;view2d, img, /cool, /colo ;, no_data_idx=where(img le 0) + ;stop + end + else : begin + print, '*** Error in save_product: not yet supported image type: ' + type + ' -> skip' + self -> print_log, 'not yet supported image type: ' + type + ' -> skip', log_mode = 32 + return + end + endcase + + ;xml_file schreiben + if not quiet then print, ' write xml file: ', xml_file + if n_elements(used_data) eq 0 then print, '*** Error in save_product: no information about used_data given' + openw, lun_1, path + xml_file, /get_lun + printf, lun_1, "<?xml version='1.0' standalone='yes'?>" + print, '<MSGProduct Name="' + name + print, '" File="' + img_file + print, '" DateOfProcessing="' + systime() + print, '" ProcessorVersion="' + version + print, '" UsedData="' + used_data + print, '" Remarks="' + remarks + '">' + printf, lun_1, '<MSGProduct Name="' + name + $ + '" File="' + img_file + $ + '" DateOfProcessing="' + systime() + $ + '" ProcessorVersion="' + version + $ + '" UsedData="' + used_data + $ + '" Remarks="' + remarks + '">' + printf, lun_1, ' <ImageInformation Type="' + img_type + $ + '" BitPerPixel="' + strcompress(string(depth, format = '(i2)'), /remove_all) + $ + '" ImgChannels="' + ((true ne 0) ? '3' : '1') + $ + '" Quality="' + strcompress(string(img_quality, format = '(f15.5)'), /remove_all) + $ + '" Scale="' + strcompress(string(scale, format = '(f15.5)'), /remove_all) + $ + '" X0="' + strcompress(string(x_0, format = '(i5)'), /remove_all) + $ + '" Y0="' + strcompress(string(y_0, format = '(i5)'), /remove_all) + $ + '" X1="' + strcompress(string(x_1, format = '(i5)'), /remove_all) + $ + '" Y1="' + strcompress(string(y_1, format = '(i5)'), /remove_all) + $ + '">' + print, ' </ImageInformation>' + print, ' <DataRepresentation Slope="' + strcompress(string(slope, format = '(f15.5)'), /remove_all) + $ + '" Offset="' + strcompress(string(offset, format = '(f15.5)'), /remove_all) + $ + '" NoDataValueRaw="' + strcompress(string(nodata_value_raw, format = '(f15.5)'), /remove_all) + $ + '" NoDataValuePhysical="' + strcompress(string(nodata_value_phy, format = '(f15.5)'), /remove_all) + $ + '" PhysicalUnit="' + unit + $ + '" Accuracy="' + strcompress(string(accuracy, format = '(f15.5)'), /remove_all) + $ + '">' + print, ' </DataRepresentation>' + print, '</MSGProduct>' + printf, lun_1, ' </DataRepresentation>' + printf, lun_1, '</MSGProduct>' + printf, lun_1, ' </ImageInformation>' + printf, lun_1, ' <DataRepresentation Slope="' + strcompress(string(slope, format = '(f15.5)'), /remove_all) + $ + '" Offset="' + strcompress(string(offset, format = '(f15.5)'), /remove_all) + $ + '" NoDataValueRaw="' + strcompress(string(nodata_value_raw, format = '(f15.5)'), /remove_all) + $ + '" NoDataValuePhysical="' + strcompress(string(nodata_value_phy, format = '(f15.5)'), /remove_all) + $ + '" PhysicalUnit="' + unit + $ + '" Accuracy="' + strcompress(string(accuracy, format = '(f15.5)'), /remove_all) + $ + '">' + printf, lun_1, ' </DataRepresentation>' + printf, lun_1, '</MSGProduct>' + free_lun, lun_1 + + endif else begin + print, '*** Error in save_product: can not write product to harddrive: ' + path + img_file + ' -> skip' + self -> print_log, '*** Error in save_product:can not write product to harddrive: ' + path + img_file + ' -> skip', log_mode = 32 + return + endelse + catch, /cancel + +end +;-------------------------------------------------------------------------------------------------------------------------------------- diff --git a/cws_read__transformit.pro b/cws_read__transformit.pro index d387e3f711daa247623ff142b683735307c97abe..1483c53193c054a08e2090e37576478459178012 100644 --- a/cws_read__transformit.pro +++ b/cws_read__transformit.pro @@ -2,7 +2,7 @@ ;+ ; procedure: ; cws_read::transformIt, force=force -; :author: Andy Walther +; :author: Andi Walther ; ; :Keywords: @@ -17,7 +17,7 @@ ;- PRO cws_read::transformIt, version=version, force=force - + version = 'ICWG2' COMMON FEEDBACK, quiet, verbose, debug verbose =1 IF verbose THEN PRINT, '... transform file format for ', self.product $ diff --git a/test_icwg2.pro b/test_icwg2.pro index c1e4aa7539f6073d10be269a5bdb65734ef4a525..6030f40daf06b1f5f23d39a9acf37fc0605bbcd3 100644 --- a/test_icwg2.pro +++ b/test_icwg2.pro @@ -2,7 +2,7 @@ pro test_icwg2 o=cws_read() -o.set_product,'cmb' +o.set_product,'ctp' o.set_group,'TPS' o.transformIt o.set_date,2008,6,13,12,0