[seek-dev] RE: garp tools

Bertram Ludaescher ludaesch at sdsc.edu
Sat Jan 31 09:57:03 PST 2004


Hi all:

Yes, integrating systems such as ArcView would be a great value-added
benefit of using Kepler. Indeed using Kepler "just" for application
integration alone is already useful.  I don't know how the specific
integration of ArcView could be done; I suspect there are both
technical and licensing issues that make this an interesting exercise
;-)

fyi: that today we have at least the following ways to integrate
"stuff" with Kepler:

-- native: write a Ptolemy actor in Java; I think there is already good
tutorial material for available from the Ptolemy group 

-- via web services: just export your application as a web service and 
provide a WSDL file -- that's it. Kepler can instantiate its generic
web service actor to your specific web service operation. Problem: web 
services don't support user interaction, GUIs etc 

-- using the browserUI actor: this is an actor which can be used both
as a display actor and as a user-input actor. Anything you can do from 
your browser (HTML incl. tables, forms etc., Javascript, SVG plug-in,
...) can be easily invoked within a Kepler pipeline. 
So the browserUI actor is a great tool for integrating existing apps
into Kepler: it can act as a data display, data entry, and
(within-the-browser) tool.

-- using the database query actor: this allows one to interface with
JDBC accessible databases. Other APIs could be added.

-- there is some support for interfacing with foreign languages
already, and more is planned. 

-- general grid service actors and EcoGrid data access actors will
provide useful as well. I think Bing will work on an EcoGrid and/or
SRB actor in the near future.

We should use the next SEEK meeting in May and possibly earlier
meetings to flesh out more of the overall Kepler extension strategy.

Bertram



>>>>> "PRS" == Pereira, Ricardo Scachetti <rpereira at ku.edu> writes:
PRS> 
PRS> 	Hi, Chad.
PRS> 	Here it is the AVENUE script that executes step (1a). It's is in
PRS> attachement.
PRS> 	Unfortunatelly, I don't know exactly how you can integrate
PRS> ArcView into Kepler. Some people at earlier SEEK-BEAM working group
PRS> meetings mentioned that ArcView would be a candidate for integration in
PRS> the pipelines, but I don't know exactly how. There is also the problem
PRS> with ArcView lincensing.
PRS> 	The best solution would be to implement the same processing step
PRS> on GRASS and have that one integrated. 
PRS> 	The best candidate for doing the job in GRASS are both r.region
PRS> and r.resample commands. Command r.region sets the extent and resolution
PRS> of the output layers and r.resample resamples the input layers. Then it
PRS> would be a matter of formating the grids to ASCII raster format.
PRS> 	I'll do some more research on the GRASS script itself.
PRS> 	Still need to come up with the code for step (1b).
PRS> 	I'll be in touch.
PRS> 	Regards,
PRS> 
PRS> Ricardo
PRS> 
PRS> 
PRS> -----Original Message-----
PRS> From: Chad Berkley [mailto:berkley at nceas.ucsb.edu] 
PRS> Sent: Monday, January 26, 2004 6:43 PM
PRS> To: Pereira, Ricardo Scachetti
PRS> Cc: seek-dev at ecoinformatics.org
PRS> Subject: Re: garp tools
PRS> 
PRS> Hi Ricardo,
PRS> 
PRS> could you send me the AVENUE script that you use for 1a?  It's possible
PRS> I can execute arcView from within kepler.  if not, i may rewrite it in
PRS> GRASS or something else.
PRS> 
PRS> for 1b, creating the actor is trivial.  would you mind sending me a c++
PRS> implementation of the vb code?  if you don't have time (I know you're
PRS> out of town), I can probably figure it out, but it would be more
PRS> efficient for you to do it.
PRS> 
PRS> thanks,
PRS> chad
PRS> 
PRS> Pereira, Ricardo Scachetti wrote:
>> Chad,
>> 
>> Here it is the expansion of "Layer Integration" step that is
PRS> SPECIFIC 
>> to GARP nowadays.
>> I'm assuming you are looking at the file 
>> /seek/projects/beam/niche-modeling/garp/GarpPipelineNative.jpg on CVS.
>> Correct me if I'm wrong.
>> In your current implementation of GARP pipeline in Kepler, you
PRS> have 
>> already implemented the analytical steps (2), (3) and (4). Now you 
>> will implement A.S. #(1): "Env Layer cliping and resampling".
>> Here is the procedure, broke down to smaller analytical steps.
>> Let me know if you need something else.
>> 
>> 
>> ==================================================================
>> 1a) Layer clipping and resampling
>> 
>> Description:
>> In this processing step, a set of environmental layers
PRS> is clipped 
>> and resampled so they all have the same extent and cell size.
>> It is assumed that all layers are stored in the same
PRS> coordinate 
>> system.
>> Currently, this step is executed inside ArcView 3.x
PRS> using an AVENUE 
>> script.
>> Can be implemented in other GIS, but that will require
PRS> development 
>> of other custom scripts that comply to GIS available API.
>> 
>> Inputs: 
>> - LAYERS: set of environmental layers (files) that will
PRS> be used in 
>> the modeling (in a format that your GIS API can read);
>> - EXT: extent of the area of interest. Defined as 2 pair
PRS> of 
>> coordinates, i.e., a (x, y) pair for upper-left corner and another for
PRS> 
>> the bottom-right corner of the area of interest;
>> - CS: cell size to be used in the analysis
>> 
>> Outputs:
>> - set of layers in ESRI ASCII Raster Grid format, all
PRS> with the same 
>> extent (EXT) and cell size (CS). Note: EXT and CS need to be exactly 
>> the same for all layers, to the last decimal digit!!
>> 
>> 
>> Note: This step can probably be broken down to smaller pieces,
PRS> but 
>> that will depend on the GIS API you are using.
>> 
>> 
>> 
>> ==================================================================
>> 1b) Scaling of environmental layer values and packaging of GARP 
>> dataset
>> 
>> Description:
>> In this step, each ASCII raster grid obtained in step 1a
PRS> is 
>> converted to native GARP dataset format. Each layer has its cell 
>> values scaled to fit a byte. This is a requirement for the original 
>> GARP algorithm as designed by David Stockwell.
>> Each layer has its maximum (max) and minimum (min) cell
PRS> value 
>> extracted. The following expression is then applied to each of the 
>> layer cells:
>> 
>> scale_coeficient = (max - min) / 253;
>> new_scaled_value = (int)( ( (original_value -
>> min) / scale_coeficient) + 1);
>> 
>> The reverse operation is (just FYI): 			
>> 
>> new_original_value = ( ( (double)
>> new_scaled_value - 1) * scale_coeficient ) + min
>> 
>> After scaling occurs, the cell values assume values
PRS> between 1 and 
>> 254 (0 and 255 are reserved values)
>> This step is currently implemented in a separate Visual
PRS> Basic 
>> application. Below there is the source code (in VB!!) that does the 
>> processing. I can translate it to C++ when you have the actor and the 
>> right method signitures (JNI) in place.
>> 
>> Input:
>> - set of ASCII Raster grids from step 1a
>> 
>> Outputs:
>> - a directory containing: i) the environmental layer
PRS> scaled values 
>> (scaled to fit a byte), one .raw file per layer;
>> - a GARP dataset file (.dxl) containing the description
PRS> (original 
>> values, extent, etc) of each layer in the dataset.
>> 
>> 
>> ==================================================================
>> 
>> 
>> 
>> Visual Basic Code (extract from VB Application source code)
>> 
>> Dim m_oDS As GARPLib.CoEnvLayerSet  ' VB version of EnvLayerSet cpp 
>> object
>> Dim m_oLyr As GARPLib.CoEnvLayer    ' VB version of EnvLayer cpp
PRS> object
>> Dim strFile as String               ' ASCII Raster file name (without
>> extention .asc)
>> 
>> Set m_oLyr = CreateObject("Garp.EnvLayer")
>> m_oLyr.strId = strFile
>> m_oLyr.strTitle = strFile
>> m_oLyr.strFilename = ""
>> m_oLyr.strMatrixFilename = strFile & ".raw"
>> Call m_oLyr.useRawByteMatrixInDisk
>> Call m_oLyr.fromAsciiGrid(strFile * ".asc")
>> 
>> If (LCase(strFile) <> "mask") Then
>> Call m_oDS.Add(m_oLyr)
>> Else
>> Call m_oDS.addAsMask(m_oLyr)
>> End If
>> 
>> ' check for error: It will raise an error if extents or cell size 
>> ' of this layer is different from the other layers
>> If m_oDS.intLastErrorCode <> 0 Then
>> Call Err.Raise(1, "Ascii Raster Grid Batch Load", _
>> "Could not add layer to the dataset """ & oFile.Name &
PRS> _
>> """ (Garp message: " & m_oDS.strLastErrorMessage &
PRS> ")")
>> End If
>> 
>> 
>> 
>> 
>> 
>> 
>> -----Original Message-----
>> From: Chad Berkley [mailto:berkley at nceas.ucsb.edu]
>> Sent: Tuesday, January 20, 2004 12:47 PM
>> To: Ricardo Scachetti Pereira
>> Cc: seek-dev at ecoinformatics.org
>> Subject: garp tools
>> 
>> Hi Ricardo (and others),
>> 
>> Now that I think that GARP is working from kepler, matt and I thought 
>> it would be a good idea to implement the tools necessary to go from 
>> raw data to the input formats that are in the example data files.  We 
>> would like to be able to implement a pipeline that goes from a digir 
>> (ecogrid) query to a garp prediction to rendered output.
>> 
>> What I need from you (and maybe deana as well) is to know exactly what
PRS> 
>> this entails and if any of it is implemented in a way that I could 
>> harness in kepler or if I need to implement it.  I believe this 
>> process is what is grouped together under "Layer Integration" in the 
>> original ppt garp pipeline that came out of BEAM. Could you possibly 
>> expand "Layer Integration" so I know exactly what the steps are within
PRS> 
>> that process?  Are any of these processes implemented in GeoTools or
PRS> GRASS?
>> 
>> thanks,
>> chad
>> 
>> 
>> --
>> -----------------------
>> Chad Berkley
>> National Center for
>> Ecological Analysis
>> and Synthesis (NCEAS)
>> berkley at nceas.ucsb.edu
>> -----------------------
>> 
PRS> 
PRS> 
PRS> --
PRS> -----------------------
PRS> Chad Berkley
PRS> National Center for
PRS> Ecological Analysis
PRS> and Synthesis (NCEAS)
PRS> berkley at nceas.ucsb.edu
PRS> -----------------------
PRS> 
PRS> 
PRS> '========================================================================================================================
PRS> ' Name    : clipresgrids.ave
PRS> ' Version : 2.1
PRS> ' Author  : Blaine Ray
PRS> ' Email   : blaine at ukans.edu
PRS> ' Date    : June 28, 2001
PRS> ' Desc.   : Clips one or more grids in the view by a selected rectangle graphic and stores the clipped grids
PRS> '           in a new directory.
PRS> ' Input   : One or more overlaying grids and one of the following: (1) a selected overlaying rectangle graphic, 
PRS> '           (
PRS> 2) a polygon shapefile, or (3) a set of user-specified coordinates.
PRS> ' Output  : Grids of the same file name that are clipped, resampled, and stored in a new directory as defined by the user.
PRS> '========================================================================================================================
PRS> 
PRS> 
PRS> ' Title for pop-up message boxes
PRS> theHeader = \"Clip and Resample Grid Themes\"
PRS> 
PRS> ' Make sure the Spatial Analyst extension is loaded
PRS> exttest=Extension.Find(\"Spatial Analyst\")
PRS> if (exttest=NIL) then
PRS>   m
PRS> sgbox.error(\"You must have the Spatial Analyst extension loaded\",theHeader)
PRS>   Exit
PRS> end
PRS> 
PRS> ' Get the path to save new grids to
PRS> newpath = MsgBox.Input(\"Input the path to which the newly clipped grids will be placed\",theHeader,\"c:\\temp\\\")
PRS> if (newpath = NIL) then
PRS>   msgbox.error(\"You must enter a valid path.\",theHeader)
PRS>   Exit
PRS> end
PRS> 
PRS> ' Get the current view
PRS> theView = av.GetActiveDoc
PRS> 
PRS> ' Get themes in the current view
PRS> theThemes = theView.getthemes
PRS> 
PRS> ' Make empty grid and feature theme lists
PRS> gridThemes = List.Make
PRS> shpThe
PRS> mes = List.Make
PRS> 
PRS> ' Populate grid and feature theme lists
PRS> For each i in theThemes
PRS>   if (i.is(GTHEME)) then
PRS>     gridThemes.add(i)
PRS>   elseif (i.is(FTHEME)) Then
PRS>     shpThemes.add(i)
PRS>   end
PRS> End
PRS> 
PRS> ' Make sure that at least one grid is present in the view
PRS> if (gridThemes.Count = 0) then
PRS>   msgbox.error(\"At least one grid theme must be present in the view.\",theHeader)
PRS>   exit
PRS> end
PRS> 
PRS> ' Select an input type to use as a clip layer
PRS> inputtype = MsgBox.ChoiceAsString({\"Selected Rectangle Graphic\",\"User-Specified Coordinates\",\"P
PRS> olygon Theme\"},\"Select a clip item.\",theheader)
PRS> if (inputtype = NIL) then
PRS>   msgbox.error(\"You must select a clipping method.\",theHeader)
PRS>   Exit
PRS> end
PRS> 
PRS> ' Initialize boolean variable maskflag as being false (no analysis mask is being used)
PRS> maskflag = FALSE
PRS> 
PRS> ' Determine the input type to initialize the clipping boundary
PRS> If (inputtype = \"Selected Rectangle Graphic\") THEN
PRS>   
PRS>   ' Get the selected rectangle
PRS>   theGraphics=theview.getgraphics
PRS>   theSelg=theGraphics.getselected
PRS>   
PRS>   ' Make sure a graphic is selected
PRS>   i
PRS> f (theSelG.count=0) then
PRS>     msgbox.error(\"Must select a single graphic in the view.\",theHeader)
PRS>     EXIT
PRS>   end
PRS>   
PRS>   ' Make sure only one graphic is selected
PRS>   if (theSelG.count>1) then
PRS>     msgbox.error(\"More than one graphic is selected, select only one\",theHeader)
PRS>     EXIT
PRS>   end
PRS>   
PRS>   ' Make sure the selected graphic is a rectangle
PRS>   agraphic=thegraphics.getselected.get(0)
PRS>   if (agraphic.getshape.is(rect).not) then
PRS>     msgbox.error(\"Shape must be a rectangle\",theHeader)
PRS>     EXIT
PRS>   end
PRS>   
PRS>   ' Set theRect as
PRS>  the rectangle shape
PRS>   theRect = agraphic.getshape
PRS> 
PRS> ElseIf (inputtype = \"User-Specified Coordinates\") THEN
PRS> 
PRS>   ' Ask user to enter bounding rectangle inputs
PRS>   llx = msgbox.input(\"Input the llx (origin) coordinate.\",theheader,\"\").asnumber
PRS>   lly = msgbox.input(\"Input the lly (origin) coordinate.\",theheader,\"\").asnumber
PRS>   urx = msgbox.input(\"Input the urx coordinate.\",theheader,\"\").asnumber      
PRS>   ury = msgbox.input(\"Input the ury coordinate.\",theheader,\"\").asnumber
PRS> 
PRS>   ' Set theRect as the four user inputs
PRS>   t
PRS> heRect = rect.makeXY(llx,lly,urx,ury)
PRS> 
PRS> ElseIf (inputtype = \"Polygon Theme\") THEN
PRS> 
PRS>   ' Make sure at least one feature theme is present in the view
PRS>   If (shpthemes.count = 0) Then
PRS>     msgbox.error(\"No shapefiles present in the view.\",theHeader)
PRS>     Exit
PRS>   End
PRS>   
PRS>   ' Ask user to select a shapefile to use as a clip layer
PRS>   cliptheme = MsgBox.list(shpthemes,\"Select a shapefile containing a single shape for clipping.\" ,theHeader)
PRS>   if (cliptheme = NIL) then
PRS>     msgbox.error(\"You must select a clip shapefile.\",the
PRS> Header)
PRS>     Exit
PRS>   end
PRS> 
PRS>   ' Get the feature table of the selected shape file
PRS>   clipftab = cliptheme.getftab
PRS> 
PRS>   ' Use the first shape in the shapefile as a clip shape
PRS>   clipshp = clipftab.returnvalue((clipftab.findfield(\"shape\")),0)
PRS>   
PRS>   ' Get the minimum bounding box for the clip shape
PRS>   clipextent = clipshp.returnextent
PRS> 
PRS>   ' Set theRect as the minimum bounding box for the clip shape
PRS>   theRect = clipextent
PRS>   
PRS>   ' Ask user if values outside of the shape but within the clip boundaries should be masked to a va
PRS> lue of NODATA
PRS>   maskflag = MsgBox.YesNo(\"Mask values outside of theme as NO DATA?\", theheader, TRUE)
PRS>   
PRS> End
PRS>   
PRS> ' Get the x,y coordinate for the lower left corner of clipping shape
PRS> rectorigin = theRect.returnorigin
PRS>  
PRS> ' Get the x,y coordinate for the upper right corner of clipping shape
PRS> rectsize = theRect.returnsize
PRS> rectUR = Point.make(rectorigin.getx+rectsize.getx, rectorigin.gety+rectsize.gety)
PRS> 
PRS> ' Ask user to select one or more grids to be clipped
PRS> Gridthemelist = MsgBox.multilist(gridthemes,\"Select one or m
PRS> ore grids to be clipped by the selected rectangle.\" ,theHeader)
PRS> if (Gridthemelist = NIL) then
PRS>   msgbox.error(\"You must select at least one grid theme.\",theHeader)
PRS>   Exit
PRS> end
PRS> 
PRS> ' Initialize the variable maxcellsz to 0 prior to the looping process
PRS> maxcellsz = 0
PRS> 
PRS> ' Initialize the string variable gmsg to be empty prior to the looping process
PRS> gmsg = \"\"
PRS> 
PRS> ' Get the number of grids to be clipped 
PRS> numgrids = gridthemelist.count
PRS> 
PRS> ' Generate a summary message for the grids to be clipped (describe number of rows, column
PRS> s, and cell sizes for each grid)
PRS> For each gthm in Gridthemelist
PRS> 
PRS>   ' Current grid name
PRS>   gname = gthm.getname
PRS> 
PRS>   ' Current grid
PRS>   ggrid = gthm.getgrid
PRS> 
PRS>   ' Current number of rows
PRS>   grows = ggrid.getnumrowsandcols.get(0).asstring
PRS>   
PRS>   ' Current number of columns
PRS>   gcols = ggrid.getnumrowsandcols.get(1).asstring
PRS>   
PRS>   ' Get current cell size
PRS>   gcellsz = ggrid.getcellsize
PRS> 
PRS>   ' Obtain the maximum grid cell size
PRS>   if (gcellsz > maxcellsz) then
PRS>     maxcellsz = gcellsz
PRS>   END
PRS> 
PRS>   ' Get new summary message text
PRS>   gtem
PRS> pmsg = (gname+\" > Dimensions: \"+grows+\" x \"+gcols+\"; Cell Size: \"+gcellsz.asstring+NL)  
PRS> 
PRS>   ' Add the text to the summary message
PRS>   gmsg = gmsg + gtempmsg
PRS> 
PRS> END
PRS> 
PRS> ' Output the summary message  
PRS> msgbox.Report(gmsg,theheader)
PRS> 
PRS> ' Ask user to declare if the grids should be resampled
PRS> resampleflag = MsgBox.YesNo(\"Resample grids to a different cell size?\"+NL+\"Note: All grids will be resampled to the same cell size.\", theheader, TRUE)
PRS> 
PRS> ' If grids are to be resampled, get the new sample size and the resampling method 
PRS>  
PRS> if (resampleflag) Then
PRS>   
PRS>   ' Initialize resampledone flag to be equal to FALSE
PRS>   resampledone = FALSE
PRS> 
PRS>   ' Until user finalizes output grid cell size, loop through the following processes 
PRS>   While (resampledone <> TRUE)
PRS>  
PRS>     ' Get the new cell size for the output grids
PRS>     newcellsz = MsgBox.Input(\"Input the new desired cell size (default is the maximum input grid cell size).\",theheader,maxcellsz.asstring).asnumber
PRS>     if (newcellsz = NIL) then
PRS>       msgbox.error(\"You must input a valid cell size.\",theH
PRS> eader)
PRS>       Exit
PRS>     end
PRS> 
PRS>     ' If grids are resampled, obtain new bounding coordinates (add a row/column to each side)
PRS>     If (resampleflag) Then
PRS>       
PRS>       'Get new coordinates (add a row/column to each side)
PRS>       neworgx = rectorigin.getx - newcellsz
PRS>       neworgy = rectorigin.gety - newcellsz
PRS>       newURx = rectUR.getx + newcellsz
PRS>       newURy = rectUR.gety + newcellsz
PRS>       
PRS>       ' Set the new bounding rectangle coordinates
PRS>       outrect = rect.makeXY(neworgx,neworgy,newURx,newURy)
PRS>     End
PRS>     
PRS>   
PRS>   ' Get output grid properties (number of rows, number of columns, and cell size) 
PRS>     newnumrows = ((outrect.getwidth) / newcellsz).round
PRS>     newnumcols = ((outrect.getheight) / newcellsz).round
PRS>     newnumcells = newnumrows * newnumcols
PRS>     
PRS>     ' Convert afforementioned values to strings
PRS>     newstrrows = newnumrows.asstring
PRS>     newstrcols = newnumcols.asstring
PRS>     newstrcells = newnumcells.asstring
PRS>     
PRS>     ' Allow user to determine if the current grid cell size is appropriate
PRS>     resampledone = MsgBox.Ye
PRS> sNo(\"Output resampled grid dimensions: \"+newstrrows+\" x \"+newstrcols+NL+\"Total number of cells: \"+newstrcells+NL+\"Is this OK?\",theheader,FALSE)
PRS>   
PRS>   End
PRS> End
PRS> 
PRS>     ' Make a resampling method string list
PRS>     resampletypelist = {\"Nearest Neighbor\",\"Bilinear Interpolation\",\"Cubic Convolution\"}
PRS>   
PRS>     ' Get the resampling method for the output grids
PRS>     resampletype = MsgBox.listasstring(resampletypelist,\"Select a resampling method.\",theheader)
PRS>     if (resampletype = NIL) then
PRS>       msgbox.error(\"You must select 
PRS> a resampling method.\",theHeader)
PRS>       Exit
PRS>     end
PRS>     
PRS>     ' Get a resampling enumeration type
PRS>     if (resampletype = \"Nearest Neighbor\") Then
PRS>       resampleenum = #GRID_RESTYPE_NEAREST
PRS>     elseif (resampletype = \"Bilinear Interpolation\") Then
PRS>       resampleenum = #GRID_RESTYPE_BILINEAR 
PRS>     elseif (resampletype = \"Cubic Convolution\") Then
PRS>       resampleenum = #GRID_RESTYPE_CUBIC
PRS>     End
PRS> 
PRS> 
PRS> ' Ask user to select an output file format
PRS> OutputType = MsgBox.ChoiceAsString({\"ARC/INFO Grid\", \"ASCII Grid\"},\"Select
PRS>  the file format for the output grids.\",theheader)
PRS> if (OutputType = NIL) then
PRS>   msgbox.error(\"You must select an output file format.\",theHeader)
PRS>   Exit
PRS> end
PRS> 
PRS> 
PRS> '============================================================================================================
PRS> 
PRS> ' Begin main looping process
PRS> For each i in Gridthemelist
PRS>   
PRS>   ' Get the curent grid
PRS>   theGrid = i.getgrid
PRS> 
PRS>   ' Get the current cell size
PRS>   cellsz=theGrid.getcellsize
PRS> 
PRS>   ' If resampleflag is TRUE then resample the grid 
PRS>   If (resampleflag) Then
PRS> 
PRS> 
PRS>     ' Reset the analysis environment
PRS>     Grid.Reset
PRS> 
PRS>     ' Resample the grid
PRS>     resgrid = thegrid.resample(newcellsz, resampleenum)
PRS> 
PRS>     ' Set analysis environment cell size to equal the newly declared cell size
PRS>     Grid.SetAnalysisCellSize(#GRID_ENVTYPE_VALUE, newcellsz)
PRS> 
PRS>   Else
PRS> 
PRS>     ' If resampleflag is FALSE then obtain new bounding coordinates (add a row/column to each side)
PRS>     neworgx = rectorigin.getx - cellsz
PRS>     neworgy = rectorigin.gety - cellsz
PRS>     newURx = rectUR.getx + cellsz
PRS>     newURy = re
PRS> ctUR.gety + cellsz
PRS> 
PRS>     ' Set the new bounding rectangle coordinates
PRS>     outrect = rect.makeXY(neworgx,neworgy,newURx,newURy)
PRS>     resgrid = thegrid
PRS>     
PRS>     ' Set the analysis cell size to be equal to the current grid's cell size
PRS>     Grid.SetAnalysisCellSize(#GRID_ENVTYPE_VALUE, cellsz)
PRS>   
PRS>   End
PRS> 
PRS>   ' Set a temporary Analysis Environment
PRS>   Grid.SetAnalysisExtent(#GRID_ENVTYPE_VALUE, outrect)
PRS> 
PRS>   ' Since shapes are already projected, set the projection as the NULL value
PRS>   theprojection=Prj.MakeNull
PRS> 
PRS>   ' If mas
PRS> kflag is TRUE then make a mask grid from the selected shapefile
PRS>   If (maskflag) Then
PRS>     clipgrid = Grid.MakeFromFTab (clipftab, theprojection, Nil, Nil)
PRS>     
PRS>     ' Set the analysis mask to be equal to clipgrid
PRS>     Grid.SetAnalysisMask(clipgrid)
PRS>   End
PRS>   
PRS>   ' Extract grid with new bounding coordinates from the resgrid 
PRS>   maskgrid = resgrid.ExtractbyRect(outrect, theProjection, FALSE)
PRS> 
PRS>   ' Get the source name of the current grid
PRS>   gridname = thegrid.getsrcname
PRS>   
PRS>   ' Determine which file format to save the ou
PRS> tput grid as
PRS>   If (OutputType = \"ARC/INFO Grid\") Then
PRS> 
PRS>     ' Assign the same name to the grid but save it in the previously specified path
PRS>     maskgrid.SaveDataSet((newpath+gridname.getname.asstring).asfilename)
PRS>     
PRS>     ' Create a new grid theme 
PRS>     gthm = GTheme.Make(maskgrid)
PRS>     
PRS>     ' Add new theme to the view
PRS>     theView.AddTheme(gthm)
PRS>   
PRS>     ' Turn new grid theme on
PRS>     gthm.setvisible(true)
PRS> 
PRS>   ElseIf (OutputType = \"ASCII Grid\") Then
PRS>     
PRS>     ' Save the grid as an ascii grid in the previously specif
PRS> ied path
PRS>     maskgrid.SaveAsAscii((newpath+gridname.getname.asstring+\".asc\").asfilename)    
PRS>   
PRS>   End
PRS> End
PRS> 



More information about the Seek-dev mailing list