SHARE
TWEET

Untitled

a guest Sep 20th, 2010 150 Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
  1. #!python
  2.  
  3. """
  4. This script contains tsx_quad pol reader
  5. and UnitTest
  6. """
  7.  
  8. from tsx_dual import tsx_dual,latlon2xy
  9. import unittest
  10. import osgeo.gdal as gdal
  11. import numpy as np
  12.  
  13. class tsx_quad(tsx_dual):
  14.    
  15.     def __init__(self,hh_pol="",vv_pol="",hv_pol="",vh_pol="",slc=""):
  16.         self.slc_ds = gdal.Open(slc)
  17.         self.meta_xml = slc
  18.         self.slc_transform = gdal.GCPsToGeoTransform(self.slc_ds.GetGCPs())
  19.         self.xsize = self.slc_ds.RasterXSize
  20.         self.ysize = self.slc_ds.RasterYSize
  21.        
  22.     def calc_raw_amp_list(self,point_list,win):
  23.         """
  24.        Calculates raw amplitude
  25.        e.g. all points in a field
  26.        """
  27.        
  28.         slc_transform = self.slc_transform
  29.        
  30.         hh_list = list()
  31.         hv_list = list()
  32.         vv_list = list()
  33.         vh_list = list()
  34.  
  35.         #DEBUG
  36.         #logfile = file("C:\\test_xy.txt","a")
  37.        
  38.         for point in point_list:
  39.             slc_x,slc_y = latlon2xy(point[0],point[1],slc_transform)
  40.  
  41.             #DEBUG
  42.             #logfile.write(str(slc_x)+","+str(slc_y)+",")
  43.             try:
  44.                 hh = self.slc_ds.GetRasterBand(1).ReadAsArray(slc_x,slc_y,win,win)
  45.                 hv = self.slc_ds.GetRasterBand(2).ReadAsArray(slc_x,slc_y,win,win)
  46.                 vh = self.slc_ds.GetRasterBand(3).ReadAsArray(slc_x,slc_y,win,win)
  47.                 vv = self.slc_ds.GetRasterBand(4).ReadAsArray(slc_x,slc_y,win,win)
  48.             except:
  49.                 hh = np.zeros((win,win))
  50.                 hv = np.zeros((win,win))
  51.                 vh = np.zeros((win,win))
  52.                 vv = np.zeros((win,win))
  53.                
  54.             hh_amp = np.abs(hh)
  55.             hv_amp = np.abs(hv)
  56.             vh_amp = np.abs(vh)
  57.             vv_amp = np.abs(vv)
  58.              
  59.             hh_list.append(np.mean(hh_amp))
  60.             hv_list.append(np.mean(hv_amp))
  61.             vh_list.append(np.mean(vh_amp))
  62.             vv_list.append(np.mean(vv_amp))
  63.  
  64.         #DEBUG
  65.         #logfile.close()
  66.        
  67.         return hh_list,hv_list,vh_list,vv_list
  68.    
  69.     def in_data(self,latlon):
  70.         slc_transform = self.slc_transform
  71.         x,y = latlon2xy(latlon[0],latlon[1],slc_transform)
  72.         return ((x > 0 and y > 0) and (x < self.xsize and y < self.ysize))
  73.    
  74.     def read_area(self,tlbr):
  75.         tllat = tlbr[0][0]
  76.         tllon = tlbr[0][1]
  77.         brlat = tlbr[1][0]
  78.         brlon = tlbr[1][1]
  79.        
  80.         xstart,ystart = latlon2xy(tllat,tllon,self.slc_transform)
  81.         xend,yend = latlon2xy(brlat,brlon,self.slc_transform)
  82.        
  83.         tlcornerx = min(xstart,xend)
  84.         tlcornery = min(ystart,yend)
  85.         xsize = abs(xend-xstart)
  86.         ysize = abs(yend-ystart)
  87.        
  88.         try:
  89.             return self.slc_ds.ReadAsArray(tlcornerx,tlcornery,
  90.                                            xsize,ysize)
  91.         except:
  92.             return None
  93.        
  94.  
  95. class TestTsxQuad(unittest.TestCase):
  96.     def setUp(self):
  97.         tsx_file = "C:\\imagery_temp\dims_op_oc_dfd2_205662493_1\\"+\
  98.         "TSX-1.SAR.L1B\\TSX1_SAR__SSC______SM_Q_DRA_20100419T201037"+\
  99.         "_20100419T201045\\TSX1_SAR__SSC______SM_Q_DRA_20100419T201037"+\
  100.         "_20100419T201045.xml"
  101.         self.test_loc = ((-34.534190,138.748411),)
  102.         self.fail_loc = ((0.0,0,0),)
  103.         self.test_area = ((-34.534190,138.748411),(-34.6,138.8))
  104.         self.tsx_data = tsx_quad(slc=tsx_file)
  105.        
  106.     def test_valid_gdal(self):
  107.         """
  108.        Valid GDAL dataset has been initialized
  109.        """
  110.         self.assertNotEqual(self.tsx_data.slc_ds,None)
  111.    
  112.     def test_band_count(self):
  113.         """
  114.        Quad pol data has 4-bands
  115.        """
  116.         self.assertEqual(self.tsx_data.slc_ds.RasterCount,4)
  117.        
  118.     def test_size(self):
  119.         """
  120.        Data size is non-zero
  121.        """
  122.         self.assertNotEqual(self.tsx_data.slc_ds.RasterXSize,0)
  123.         self.assertNotEqual(self.tsx_data.slc_ds.RasterYSize,0)
  124.        
  125.     def test_read(self):
  126.         """
  127.        Data can be read from GDAL
  128.        """
  129.         data = self.tsx_data.slc_ds.ReadAsArray(0,0,1,1)
  130.         self.assertEqual(data.shape[0],4)
  131.        
  132.     def test_read_outside(self):
  133.         """
  134.        Read outside extents is 0
  135.        """
  136.         mult_data = self.tsx_data.calc_raw_amp_list(self.fail_loc, 1)
  137.         self.assertNotEqual(mult_data,None)
  138.         self.assertEqual(len(mult_data),4)
  139.         self.assertEqual(len(mult_data[0]),1)
  140.         for channel in mult_data:
  141.             self.assertEqual(channel[0],0)
  142.            
  143.     def test_read_multiple(self):
  144.         """
  145.        Read in extents is non-zero [town]
  146.        """
  147.         mult_data = self.tsx_data.calc_raw_amp_list(self.test_loc, 1)
  148.         self.assertNotEqual(mult_data,None)
  149.         self.assertEqual(len(mult_data),4)
  150.         self.assertEqual(len(mult_data[0]),1)
  151.         for channel in mult_data:
  152.             self.assertNotEqual(channel[0],0)
  153.        
  154.     def test_read_area(self):
  155.         """
  156.        Read a set of pixels given topleft and bottom right
  157.        """
  158.         mult_data = self.tsx_data.read_area(self.test_area)
  159.         self.assertNotEqual(mult_data,None)
  160.         self.assertNotEqual(np.mean(mult_data),0)
  161.        
  162.     def test_in_data(self):
  163.         """
  164.        Do latlon in data bounds checking
  165.        """
  166.         self.assertTrue(self.tsx_data.in_data(self.test_loc[0]))
  167.         self.assertFalse(self.tsx_data.in_data(self.fail_loc[0]))
  168.    
  169.    
  170. if __name__=="__main__":
  171.     unittest.main()
RAW Paste Data
Top