Dear Ansley,Thanks for your code. I have one query which I forgot to mention in my earlier post that in the attached new code (temp.jnl) if I want to take average of xgp and ygp so that it will list k=1;105 values as the output.What I tried is:
let cne_data_2000 = if gp eq 1 then gp*Optical_Depth_Land_And_Ocean_Mean[d=1,x=65:95@ave, y=5:40@ave,k=1:105] but not working!any change? Please suggest...On Tue, Jan 13, 2015 at 4:02 AM, Ansley Manke <ansley.b.manke@xxxxxxxx> wrote:Hi,
The pt_in_poly function makes a variable on the same XY grid as the incoming variable, so it's just marking whether grid points lie inside or outside the polygon.
Here are the commands in your script,
let xgp = {77,78,78,81,81,82,85,85,88,88,86,86,87,87,86,86,82,82,83,83,84,84,83,83,79,79,78,77,77}First, I would suggest that you make your polygon a closed polygon by repeating the first point at the end of the lists. This ensures that the function will work correctly. It works by drawing lines across the polygon and seeing how many times the line crosses the outline of the polygon. It might work fine with the unclosed polygon, but this will make it certain, and it means that you can draw the polygon on a map with "plot/vs/line xgp, ygp" .
let ygp = {31,31,30,30,29,29,28,28,27,24,24,24,23,20,20,19,19,20,20,21,21,24,24,25,25,27,27,27,27}
let gp = pt_in_poly(rf[d=1],xgp, ygp)
let cne_data_2000 = if gp eq 1 then gp*rf[d=1,l=18050:18171@ave]
let xgp = {77,78,78,81,81,82,85,85,88,88,86,86,87,87,86,86,82,82,83,83,84,84,83,83,79,79,78,77,77,77}Now, I assume that what you want is only the data inside the polygon, listed with their longitude and latitude values. You'll need to make the result of the function into a 1-D list, define variables with the corresponding longitude and latitude data, and then list out just the "lon, lat, data" for the points within the polygon.
let ygp = {31,31,30,30,29,29,28,28,27,24,24,24,23,20,20,19,19,20,20,21,21,24,24,25,25,27,27,27,27,31}
! Define 2-D arrays, with the longitude values and the latitude
! values at each location of the grid.
! (Try a SHADE plot of these variables to see what we are doing here)
let xpts = x[gx=rf[d=1]] + 0*y[gy=rf[d=1]]
let ypts = 0*x[gx=rf[d=1]] + y[gy=rf[d=1]]! 1-Dimensional lists with the locations and all of the data in the
! variable cne_data_2000
let xlist = xsequence(xpts)
let ylist = xsequence(ypts)
let cne_data_2000_list = xsequence(cne_data_2000)
! We want just the locations where cne_data_2000 is valid.
! To pick out just some data in a variable use a mask.
let mask = if cne_data_2000_list then 1
let cne_out = mask * cne_data_2000_list
let x_out = mask * xlist
let y_out = mask * ylist
! The COMPRESS functions will move the valid values of these
! to the start of the list. Then we can list just the good data
! at the start of the 1-D lists.
! Make Ferret definitions for these variables, so the listing will have
! a useful header. You will want to set up the right units and a
! descriptive title for the variable.
let/title=longitude/units=degrees_east lon = compressi(x_out)
let/title=latitude/units=degrees_north lat = compressi(y_out)
let/units=" "/title=" " rf_in_region = compressi(cne_out)
let npts = mask[i=@ngd]
list/i=1:`npts` lon, lat, rf_in_region
On 1/11/2015 8:22 PM, Nitin Patil wrote:
Dear ferret users,
In the attached script I want to extract cne_data_2000 into ascii but what is happening is when I am doing:
list cne_data_2000 it is writing all the grid points in the domain. I require the points/values to be listed xgp and ygp only. Anyone knows where I am making mistake.
--Regards,Nitin Patil