Catalog

1. Class content

1. obtain Modis Swath data ( Only get here Lat、Lon、Aod( Aerosol thickness ) Three data sets and aod Two properties of data set ), Also on aod Simple processing of data

2. Call the secondary development interface to Lat、Lon Data is created as raw materials Glt file

3. Call the secondary development interface to Glt File as reference , Yes Aod Data is re projected and output

2. ENVI Practice and IDL Programming

1. Class content

Batch reprojection

Here is the overall idea of re projection ：

1. obtain Modis Swath data ( Only get here Lat、Lon、Aod( Aerosol thickness ) Three data sets and aod Two properties of data set ), Also on aod Simple processing of data

2. Call the secondary development interface to Lat、Lon Data is created as raw materials Glt file

3. Call the secondary development interface to Glt File as reference , Yes Aod Data is re projected and output

The main reason here is that the secondary development interface is relatively strange , In fact, it's the calling package. I feel .

In addition, it is relatively strange to use the secondary interface ENVI First perform the operation of re projection, and then ENVI Re projection operation and IDL Code for comparison , You will find code and ENVI The re projection operation of has many similarities , This will help you understand the operation and operation of the code

2. ENVI Practice and IDL Programming

function get_hdf4_ds, file_path, ds_name
  ;  Building functions ,  Used to obtain hdf4 Data of data set in file 
  
  ;  get files id
  file_id = hdf_sd_start(file_path, /read)
  
  ;  Get data set index
  ds_index = hdf_sd_nametoindex(file_id, ds_name)
  
  ;  Get the id
  ds_id = hdf_sd_select(file_id, ds_index)  ;  Incoming data set id And the name of the dataset 
  
  ;  Get the data of the dataset 
  hdf_sd_getdata, ds_id, ds_data
  ;  Incoming data set id, Pass in a variable to receive the data of the data set obtained by the function 
  
  ;  Close data sets and files that have been opened 
  hdf_sd_endaccess, ds_id
  hdf_sd_end, file_id
  
  ;  Return value 
  return, ds_data
end


function get_hdf4_att_info, file_path, ds_path, att_name
  ;  Building functions :  Used to get hdf4 Relevant attributes of the data set of the file 
  
  ;  get files id
  file_id = hdf_sd_start(file_path, /read)
  
  ;  Get the index
  ds_index = hdf_sd_nametoindex(file_id, ds_path)
  
  ;  Get the id
  ds_id = hdf_sd_select(file_id, ds_index)
  
  ;  Get the properties under the dataset index
  att_index = hdf_sd_attrfind(ds_id, att_name)
  
  ;  Get some attribute content of the dataset 
  hdf_sd_attrinfo, ds_id, att_index, data=att_data
  
  ;  Close datasets and files 
  hdf_sd_endaccess, ds_id
  hdf_sd_end, file_id
  
  ;  Return value 
  return, att_data
end


pro week_four_study1
  ;  Yes modis product (HDF4) Re project ( Two steps ,  establish glt file , borrow glt File projection )
  
  ;  Record the start time 
  start = systime(1)
  
  ;  Because we need to use envi Of api Interface , So you need to declare 
  compile_opt idl2
  envi,/restore_base_save_files
  envi_batch_init
  
  ;  First of all get modis In the document lon、lat( The first two data are used to create glt That is, geographic information lookup table )、aod data , And stored in the hard disk , Convenient for follow-up envi Of api To call the file 
  
  ;  route 
  in_path = 'D:/IDL_program/experiment_data/chapter_3/modis_swath'
  out_path = 'D:/IDL_program/experiment_data/chapter_3/modis_swath/geo_out'
  
  ;  If out_path directory does not exist , So create 
  if file_test(out_path, /directory) eq 0 then begin
    file_mkdir, out_path
  endif
  
  ;  obtain in_path All in modis data ( That is to get hdf4 file )
  ;  Get all hdf4 Path to file ( With array returns )
  file_path_array = file_search(in_path, '*.hdf')  ;  Incoming directory , Pass in what you need to find hdf4 The suffix of the file is .hdf
  ;  Get all under this directory hdf4 Number of files 
  file_count = n_elements(file_path_array)
  
  
  ;  obtain lon、lat、aod Data sets and aod Several properties of data sets 
  for file_i = 0, file_count - 1 do begin
    ;  The path of the file under this cycle 
    file_path = file_path_array[file_i]
    ;  Get dataset data 
    lon_data = get_hdf4_ds(file_path, 'Longitude')
    lat_data = get_hdf4_ds(file_path, 'Latitude')
    aod_data = get_hdf4_ds(file_path, 'Image_Optical_Depth_Land_And_Ocean')
    
    ;  obtain aod Properties of the dataset (_FillValue, scale_factor)
    aod_fv = get_hdf4_att_info(file_path, 'Image_Optical_Depth_Land_And_Ocean', '_FillValue')
    aod_sf = get_hdf4_att_info(file_path, 'Image_Optical_Depth_Land_And_Ocean', 'scale_factor')
    ;  We need to pay attention to , The attribute value returned here is in the form of array array( Even if it has only one float Count ), Therefore, attention should be paid to the following use 
    
    ;  Yes aod Data processing 
    aod_data = (aod_data ne aod_fv[0]) * aod_fv[0] * aod_sf[0]  ;  Be careful aod_fv and aod_sf Are arrays 
    
    ;  The output path ( Maybe you'll find out , Each cycle will overwrite the previous directory , But don't worry , We will use the data stored until now in this cycle )
    lon_path = out_path + '/lon_out.tiff'
    lat_path = out_path + '/lat_out.tiff'
    aod_path = out_path + '/aod_out.tiff'
    
    ;  Store data from three data sets 
    write_tiff, lon_path, lon_data, /float
    write_tiff, lat_path, lat_data, /float
    write_tiff, aod_path, aod_data, /float
    
    ;  Give the stored file to envi Of api call ( Or we tell api Where are the documents he needs )
    envi_open_file, lon_path, r_fid=lon_id
    envi_open_file, lat_path, r_fid=lat_id
    envi_open_file, aod_path, r_fid=aod_id
    ;  The path to the incoming storage data set , Pass parameters by keywords (r_fid=) Get the data set ( Or stored tiff file ) Of id( With variable lon_id Waiting to receive )
    
    ;  Output glt Path to file 
    out_path_glt = out_path + '/glt.img'
    out_path_glt_hdr = out_path + '/glt.hdr'  ;  Is it OK without this , Try it later 
    
    ;  Input and output projection information 
    in_prj = envi_proj_create(/GEOGRAPHIC)  ;  The input is latitude and longitude 
    out_prj = envi_proj_create(/GEOGRAPHIC)  ;  The output is also longitude and latitude 
    
    ;  establish glt file 
    envi_glt_doit, $  ; $  Means line break 
      i_proj=in_prj, x_fid=lon_id, y_fid=lat_id, x_pos=0, y_pos=0, $  ;  establish glt Relevant information that needs to be passed in 
      o_proj=out_prj, pixel_size=pixel_size, rotation=0.0, $  ;  Output glt Relevant information that needs to be passed in 
      out_name=out_path_glt, r_fid=glt_id  ;  This is also , The output path and the returned path are created glt Of documents id
      
    ;  Output the path of the re projection result 
    out_georef_path = out_path + '/' + file_basename(file_path, '.hdf') + '_georef.img'
    out_georef_hdr_path = out_path + '/' + file_basename(file_path, '.hdf') + '_georef.hdr'
    
    ;  Use glt Re project the file 
    envi_georef_from_glt_doit, $
      glt_fid=glt_id, pos=0, $  ;  Pass in glt Of documents id And the number of layers ( Because there is only one layer to fill 0 that will do )
      fid=aod_id, $  ;  Pass in the to be projected aod Of documents id
      out_name=out_georef_path, r_fid=georef_id  ;  Pass in the path where the projected file is stored , And return the file id
      
    ;  Close open files 
    envi_file_mng, id=lon_id, /remove
    envi_file_mng, id=lat_id, /remove
    envi_file_mng, id=aod_id, /remove
    envi_file_mng, id=glt_id, /remove
    envi_file_mng, id=georef_id, /remove
    
    ;  Some files need to be deleted in the middle of the re projection , Because we only need to re project the file , As for the latitude and longitude file 、glt If you want to keep the file, please don't delete 
    file_delete, [lon_path, lat_path, aod_path, out_path_glt, out_path_glt_hdr]
    
    ;  Record the end time 
    stop = systime(1)
  endfor
  print, stop - start
end

Compile operation ：