correct the color bar with the color of countor line

46 views
Skip to first unread message

Khaled Tahoon

unread,
Nov 21, 2025, 5:40:35 PMNov 21
to idl-pvwave
Hello,
I have an irregular data of lon, lat and values parameters say  r_mean over the middle east. How you can plot them as a counter lines over Egypt only with a color bar represent only the contour line over Egypt
in other words
why my color bar is showing the full Middle-East range and how to restrict it to Egypt only
So, I wort the following code
;==========================
pro contour_plot, postscript=postscript
file='C:\Users\.....\data.txt'

nrec_m = file_lines(file)
data_m=fltarr(3,nrec_m)
openr,lun,file,/get_lun
readf,lun,data_m
free_lun,lun
Lon       = reform(data_m[0,*])
Lat       = reform(data_m[1,*])
r_mean    = reform(data_m[2,*])

===================
lat_g= Scale_Vector(Findgen(100), 21.99,31.67)     ;          give an array of 100 element
lon_g = Scale_Vector(Findgen(100), 24.69,39.2);36.91)   ;          give an array of 100 element
r_pw = Scale_Vector(Findgen(100), min(r_mean), max(r_mean))    ;          give an array of 100  element

limit=[min(lat_g), min(lon_g), max(lat_g), max(lon_g)]
latsubs  = Value_Locate(lat_g, [limit[0],limit[2]])  ;   Idl print 0     100
lonsubs = Value_Locate(lon_g, [limit[1],limit[3]]) ;   Idl print 0     100

lat_sub  = lat_g[latsubs[0]:latsubs[1]]       ; Idl print the lat interpulate over the rang from 0 to 100 and become 100 values
lon_sub  = lon_g[lonsubs[0]:lonsubs[1]]    ; Idl print the lat interpulate over the rang from 0 to 100 and become 100 values

;============ convert irregular data set to a grid data

GRID_INPUT, lon, lat, r_mean , xSorted, ySorted, dataSorted

; Initialize the grid parameters.
gridSize = [100,100]
; Use the equation of a straight line and the grid parameters to
; determine the x of the resulting grid.
       slope = (MAX(xSorted) - MIN(xSorted))/(gridSize[0] - 1)
       intercept = MIN(xSorted)
       xGrid = (slope*FINDGEN(gridSize[0])) + intercept

; Use the equation of a straight line and the grid parameters to
; determine the y of the resulting grid.
    slope = (MAX(ySorted) - MIN(ySorted))/(gridSize[1] - 1)
    intercept = MIN(ySorted)
    yGrid = (slope*FINDGEN(gridSize[1])) + intercept

; Grid the data with the Radial Basis Function method.
   grid = GRIDDATA(xSorted, ySorted, dataSorted, $
   DIMENSION = gridSize, METHOD = 'RadialBasisFunction')

 position=[0.1,0.1,0.9,0.84]
; Open a display sumdow to  make a mask.
   cgDisplay, 800, 675, Color='black',  /Free, /Pixmap
 ; Set up the map projection of the US continent.
                                                 ;cgLoadCT, 33, NColors=10, Bottom=1
   cgLoadCT, 27, NColors=15, Bottom=1              ;/Reverse, /Brewer, NColors=15
   TVLCT, cgColor('white', /Triple), 254
    nLevels= 15
    c_colors = Indgen(nLevels)+1

     inc=(Max(grid) - Min(grid)) / nLevels
     Levels = Min(grid) + inc * FINDGEN(ROUND((Max(grid) - Min(grid)) / inc) + 1)

    tickname=String(nLevels)


   cgMap_Set,map_center_lat, map_center_lon, 0,  LIMIT=[21.0,24.0,32.0,37.5],   $
      /MERCATOR, /NOERASE, POSITION=position

 ; Draw the Egypt outline in white on the black background. Make a mask.
   shapefile = 'C:\Users\.........\egy_admbnda_adm0_capmas_itos_20170421.shp'

cgDrawShapes, shapefile, AttrName='ADM0_EN',  Thick=1, $
       Fill = 1, Colors='white', MapCoord=map, /Projected_XY ; if Colors='black' then the nile river will appear ,,,,  FColors=['white'],
   mask = TVRD() GT 0 ; white

cgDisplay, 800, 675, PIXMAP=Keyword_Set(postScript), /FREE

      ;**************ncolors =n_elements(contourLevels); 20
  cgLoadCT, 27, /Brewer, /Reverse, NColors=15, Bottom=1

 cgContour, grid, xgrid, ygrid, /OVERPLOT, LEVELS=Levels, Position=position, $
        label=2,/FILL, /Outline, COLOR='black',c_color=c_colors,c_charsize=1.5
;==========================
 ; Take a snapshot of what is currently in the sumdow.
   contourImage = cgSnapshot(TRUE=3)

   ; Multiply the contourImage times the mask.
   zeros = Where(mask EQ 0)
   background_color = cgColor('white', /Triple)
   FOR j=0,2 DO BEGIN
        r_pw  = contourImage[*,*,j] * mask
        r_pw[zeros] = background_color[j]
        contourImage[*,*,j] = r_pw
   ENDFOR


;;================
; Set up the PostScript device, if needed.
   IF Keyword_Set(postscript) THEN cgPS_Open,FILENAME='C:\Users\....\out_vis.eps'

cgImage, contourImage


 ; Redisplay the image.
cgtext,10500,2550, /device, '!7Sudan!X', CHARSIZE=1.5,ORIENTATION=0.0 ; , 350,80,
cgText,3250,10000,/device,  '!7Libya!X', CHARSIZE=1.5, ORIENTATION=-90.0
cgText,15750,10000,/device, '!7Red Sea!X', CHARSIZE=1.5,ORIENTATION=-60.0
cgText,10500,15750,/device, '!7Mediterranean Sea!X', CHARSIZE=1.5,ORIENTATION=0.0
cgtext,10500,1000, /device, '!7Longitude!X', CHARSIZE=1.5,ORIENTATION=0.0
cgText,12550,10000, /device, '!7Latitude!X', CHARSIZE=1.5, ORIENTATION=-90.0

   xstart =0.9
   xend = 0.95
   yend = 0.15
   ystart =0.95
c_colors = Indgen(nLevels)+1
   cgColorbar, Divisions=c_colors, Range=[Min(r_mean), Max(r_mean)], /Discrete, Format='(F4.1)', $
    Position=[0.1, 0.90, 0.9, 0.94], NColors=nLevels,Bottom=1, Title='PW_climate', Font=0, CHARSIZE=1,Ticknames=tickvalues, TLocation='Top'

   cgMap_Grid, LatDel = 1.0, LonDel =1.0, /NO_GRID , /Box_Axes, Color='white', LATLABEL=0;, $

   ; Close the PostScript device and create the PNG file with ImageMagick
   IF Keyword_Set(postscript) THEN cgPS_Close;, /png
   END
;==========================
I got the followVis_data.jpging map
data.txt
egy_admbnda_adm0_capmas_itos_20170421.shp
egy_admbnda_adm0_capmas_itos_20170421.dbf
egy_admbnda_adm0_capmas_itos_20170421.shx
egy_admbnda_adm0_capmas_itos_20170421.shp.xml

markus.sc...@gmail.com

unread,
Nov 24, 2025, 12:08:21 PM (14 days ago) Nov 24
to idl-pvwave
Hi Khaled, 
your problem is not in the colorbar, but in the contour levels. The following lines are causing your issue.
    inc=(Max(grid) - Min(grid)) / nLevels
    Levels = Min(grid) + inc * FINDGEN(ROUND((Max(grid) - Min(grid)) / inc) + 1)  
Replace Max(grid) and Min(grid) with the proper limits, and it should work. You can get it by 
    max_egypt=max(grid[where(mask ne 0)], min=min_egypt)
You may need to change the order of evaluation to get mask early enough. Alternatively, just set it by hand.
I hope that helps, Markus


Reply all
Reply to author
Forward
0 new messages