>>> from liblas import srs
  >>> from liblas import point
  >>> from liblas import header

  >>> import liblas
  >>> s = srs.SRS()
  >>> s.proj4
  ''
  >>> s.set_proj4(b'+proj=utm +zone=15 +datum=WGS84 +units=m +no_defs')
  True
  >>> print s.proj4
  +proj=utm +zone=15 +datum=WGS84 +units=m +no_defs 
  >>> s.proj4 == '+proj=utm +zone=15 +datum=WGS84 +units=m +no_defs '
  True

  >>> s = srs.SRS()
  >>> s.set_userinput(b'EPSG:4326')
  True
  >>> s.proj4 == '+proj=longlat +datum=WGS84 +no_defs '
  True

  >>> from liblas import file
  >>> f = file.File('../test/data/1.2_3.las',mode='r')
  >>> s = f.header.srs
  >>> s.wkt == b"""PROJCS["NAD83 / UTM zone 15N",GEOGCS["NAD83",DATUM["North_American_Datum_1983",SPHEROID["GRS 1980",6378137,298.2572221010002,AUTHORITY["EPSG","7019"]],AUTHORITY["EPSG","6269"]],PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433],AUTHORITY["EPSG","4269"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",-93],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AUTHORITY["EPSG","26915"]]"""
  True

  >>> s2 = srs.SRS()
  >>> s2.wkt = b"""GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.01745329251994328,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4326"]]"""
  >>> p = f.read(0)
  >>> p.x
  470692.44
  >>> f.set_srs(s2)
  True
  >>> p = f.read(0)
  >>> s2.vlr_count()
  4
  >>> s2.GetVLR(0).recordlength
  64
  >>> int(round(p.x))
  -93
  >>> int(round(p.y))
  42

  >>> del f


# -93.3515625902 41.5771483954

  >>> def new_offset(old_scale, new_scale, old_offset, x):
  ...       return (new_scale*(x - old_offset) - old_scale*x)/(-1.0*old_scale)

  >>> utm_wkt = b"""PROJCS["NAD83 / UTM zone 15N",GEOGCS["NAD83",DATUM["North_American_Datum_1983",SPHEROID["GRS 1980",6378137,298.2572221010002,AUTHORITY["EPSG","7019"]],AUTHORITY["EPSG","6269"]],PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433],AUTHORITY["EPSG","4269"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",-93],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AUTHORITY["EPSG","26915"]]"""

  >>> dd_wkt = b"""GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.01745329251994328,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4326"]]"""
  >>> s_dd = srs.SRS()
  >>> s_dd.wkt = dd_wkt
  >>> s_utm = srs.SRS()
  >>> s_utm.wkt = utm_wkt

  >>> f = file.File('../test/data/1.2_3.las',mode='r')
  >>> dd_header = liblas.header.Header(handle=f.header.handle,copy=True)
  >>> utm_header = liblas.header.Header(handle=f.header.handle,copy=True)
  >>> del f

  >>> utm_header.offset
  [0.0, 0.0, 0.0]
  >>> utm_header.scale
  [0.01, 0.01, 0.01]

#  >>> new_offset(0.01, 0.0001, 0.0, 470692.44)
#  >>> utm_header.offset = [offset+1.0/0.000001 for offset in utm_header.offset]
#  >>> utm_header.offset
  [1000000.0, 1000000.0, 1000000.0]
#  >>> utm_header.scale = [0.000001,0.000001,0.000001]
  >>> utm_header.srs = s_utm

#  >>> dd_header.scale = [0.000001,0.000001,0.000001]
  >>> dd_header.srs = s_dd

  >>> f = file.File('../test/data/1.2_3.las',mode='r', header = utm_header)
  >>> int(f.header.data_offset)
  438
  >>> f.header.scale
  [0.01, 0.01, 0.01]
  >>> p = f.read(0)
  >>> origx, origy = p.x, p.y
  >>> "{:.2f}, {:.2f}".format(origx, origy)
  '470692.44, 4602888.90'
  >>> f.set_srs(s_dd)
  True
  >>> p = f.read(0)

We only get truncated values because our header scale
values are 0.01

  >>> "{:.2f}, {:.2f}".format(p.x, p.y)
  '-93.35, 41.58'

#real values
#  (-93.351562590199833, 41.577148395415108)



  >>> f_project = file.File('junk_srs_project.las',mode='w',header=dd_header)

  >>> p.header = dd_header
  >>> "{:.2f}, {:.2f}".format(p.x, p.y)
  '-93.35, 41.58'

  >>> dd_header.srs.proj4
  '+proj=longlat +datum=WGS84 +no_defs '
  >>> f_project.write(p);f_project.write(p);f_project.write(p)
  >>> f_project.close()
  >>> del f_project
  >>> f3 = file.File('junk_srs_project.las')
  >>> int(f3.header.data_offset)
  789

  >>> s_utm = srs.SRS()
  >>> s_utm.wkt = utm_wkt
  >>> p3 = f3.read(1)
  >>> int(round(p3.x)), int(round(p3.y))
  (-93, 42)
  >>> p3 = f3.read(0)
  >>> int(round(p3.x)), int(round(p3.y))
  (-93, 42)

  >>> import os
  >>> os.remove('junk_srs_project.las')


  >>> f = file.File('../test/data/srs_vertcs.las',mode='r')
  >>> s = f.header.srs
  >>> str(s.get_wkt_compoundok().decode())
  'COMPD_CS["unknown",PROJCS["WGS 84 / UTM zone 17N",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433],AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",-81],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AUTHORITY["EPSG","32617"]],VERT_CS["NAVD88 height",VERT_DATUM["North American Vertical Datum 1988",2005,AUTHORITY["EPSG","5103"],EXTENSION["PROJ4_GRIDS","g2012a_conus.gtx,g2012a_alaska.gtx,g2012a_guam.gtx,g2012a_hawaii.gtx,g2012a_puertorico.gtx,g2012a_samoa.gtx"]],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AXIS["Up",UP],AUTHORITY["EPSG","5703"]]]'

  >>> s2 = srs.SRS()
  >>> s2.wkt = b"""PROJCS["WGS 84 / UTM zone 17N",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433],AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",-81],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AUTHORITY["EPSG","32617"]]"""
  >>> s2.set_verticalcs( 5703, b'abc', 5103, 9001 )
  True
  >>> str(s2.get_wkt_compoundok().decode())
  'COMPD_CS["unknown",PROJCS["WGS 84 / UTM zone 17N",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433],AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",-81],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AUTHORITY["EPSG","32617"]],VERT_CS["NAVD88 height",VERT_DATUM["North American Vertical Datum 1988",2005,AUTHORITY["EPSG","5103"],EXTENSION["PROJ4_GRIDS","g2012a_conus.gtx,g2012a_alaska.gtx,g2012a_guam.gtx,g2012a_hawaii.gtx,g2012a_puertorico.gtx,g2012a_samoa.gtx"]],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AXIS["Up",UP],AUTHORITY["EPSG","5703"]]]'