How to plot a contour map in UTM projection?

I am trying to plot data over a regular grid(GRID_VEL50 corresponding to Northings and Eastings in meters-UTM map projection) in the form of filled contours. I want the corresponding map overlayed along with it. Following are the variables used in the code: GRID_VEL50 (dim:156X150),GRID_X (Eastings,dim:150) and GRID_Y (Northings,dim:156). Here is the script, I tried for plotting the same:

    axesm('tranmerc','MapLatLimit',[7 37.5],'MapLonLimit',[66.667 98.6666717529297]); 
    framem; gridm; axis off; tightmap;
    x11 =245895.96;  %meter, the North-West corner(66.667 deg in long)
    y11 =774370.68;  %meter, the North-West corner(7 deg in lat) 
    dx1 = 25000;     %in meter
    dy1 =-25000;     %in meter
    R=makerefmat(x11, y11, dx1, dy1);
    contourfm(GRID_VEL50,R,'r','LineWidth',0.0001)
    coast = load('coast');
    geoshow(coast.lat, coast.long, 'Color', 'black')
    contourcbar

On running this script, it gives an error, stating :"Error using constructGeoRasterReference (line 45)." The referencing array comes out to be a 3X2 array as mentioned in help documentation.

R =

                           0                    -25000
                       25000                         0
                   220895.96                 799370.68

I am unable to debug this. Is there any other method to plot contours in a specific map projection?

[Walter Roberson]

Please show the complete error message. For example is it complaining about incompatible dimensions? Is it complaining about index out of range?

[Tanvi Gupta]

Error using constructGeoRasterReference (line 45)
In combination with the number of rows in the raster grid, 156, the input referencing vector or matrix implies latitude limits that extend outside the interval [-90 90] degrees.

Error in refmatToGeoRasterReference (line 159)
R = constructGeoRasterReference(rasterSize, rasterInterpretation, ...

Error in internal.map.convertToGeoRasterRef (line 26)
R = refmatToGeoRasterReference(R, rasterSize, func_name, var_name_R, arg_pos);

Error in parseContourInputs (line 27)
R = internal.map.convertToGeoRasterRef( ...

Error in contourm (line 97)
[Z, R, levelList, pvpairs] = parseContourInputs(varargin);

Error in contourfm (line 36)
          contourm(varargin{:},'Fill','on','DefaultLineColor','black');

[Walter Roberson]

You indicate lat 7 long 66.667 as your North-West corner, but your axesm latitudes are from 7 to 37.5 which implies that lat 7 is to be a South corner rather than a North corner. Please re-check which points are the North or South extremes and please check the sign of your dx1 and dy1 .

When I calculate with the dx1 and dy1 values you give, I do not as yet see any reason why the latitude would exceed +90 to -90, but I do see that it would be outside your map axes which suggests you have incorrect bounds or incorrect signs.

[Tanvi Gupta]

Rightly pointed. I changed y11 and gave the run. Still, same error message is coming.

y11 =4152898;  % the upper left corner(37.5 deg in lat)

I am converting lat(37.5 N) and long(66.667) to utm coordinate points to get x11 and y11.

Answer 1

[Walter Roberson]

In the documentation it is not clear to me how makerefmat distinguishes between the x/y case and the long/lat case. I speculate that if you pass integral numbers that it might interpret them as x/y, so try with

      R=makerefmat(round(x11), round(y11), round(dx1), round(dy1));

[Tanvi Gupta]

Tried with integral numbers. It still shows the same error. I am using the example given on following link for making referencing array: http://in.mathworks.com/help/map/ref/makerefmat.html

If I use a referencing vector [cells/deg,North lat,West lon], e.g.[4,37.5,66.667] contours with map overlayed are generated. But I know this is not correct because it gives a wrong spread of variable over the map.