; interrogates .ref file and loads it
; called with the ModelChosen, and the ModelLongN and LatN that the .ref ought to conform to
; returns MapIDLReg, RegSize, RegName, RegColour, RegN

pro LoadRef, ModelChosen, ModelLongN, ModelLatN, MapIDLReg, RegSize, RegName, RegColour, RegN,CallFile=CallFile

ExitStatus = 0
RefFilePath = "null" & 	RefTitle = ""
CheckLongN = 0 & 	CheckLatN = 0 &		RegN = 0L
ExeCoord = 0.0D &	WyeCoord = 0.0D
MissVal = -999.0
DatumCha = ""
Int0 = 0L &	Int1 = 0L &	Int2 = 0L &	Int3 = 0L
Int4 = 0L &	Int5 = 0L &	Int6 = 0L &	Int7 = 0L & Int8 = 0L & Int9 = 0L

while (ExitStatus NE 1) do begin				; initial interrogation
 if (keyword_set(CallFile) AND RefFilePath EQ "null") then begin
  OpenStatus = -1 & RefFilePath = CallFile
  openr, lunRef, RefFilePath, /get_lun, error=OpenStatus
  if (OpenStatus NE 0) then print, "  > File cannot be opened. Enter region .ref filepath."
 endif else begin
  print, "  > Enter the filepath of the .ref file: "
  OpenStatus = -1
  while (OpenStatus NE 0) do begin
   read, RefFilePath
   RefFilePath = strtrim(RefFilePath,2)
   openr, lunRef, RefFilePath, /get_lun, error=OpenStatus
   if (OpenStatus NE 0) then print, "  > File cannot be opened. Re-enter filepath."
  endwhile  
 endelse
 
 readf, lunRef, RefTitle,  			format='(a40)'			; title
 readf, lunRef, RegN, CheckLongN, CheckLatN, 	format='(I6,2(I4))'		; dims

 if (CheckLongN EQ ModelLongN AND CheckLatN EQ ModelLatN) then begin
   ExitStatus = 1
 endif else begin
   print, "  > File is on the wrong model grid. Re-enter filepath."
   free_lun, lunRef
 endelse
endwhile

print, "  > ", strtrim(RefTitle,2)				; print to show success in identifying file

MapIDLReg       = lonarr (ModelLongN, ModelLatN)
MapIDLReg [*,*] = long(MissVal)
RegSize         = lonarr (RegN)
RegSize [*]     = 0L
RegName         = strarr (RegN)
RegName [*]     = ""
RegColour       = intarr (RegN)
RegColour [*]   = 255

for XReg = 0L, (RegN-1L) do begin					; load RegName
    readf, lunRef, DatumCha, format="(A20)"
    RegName (XReg) = DatumCha
endfor

if ((float(ModelLongN)/8.0) EQ fix(float(ModelLongN)/8.0)) then begin
 LongEightN = round (ModelLongN / 8)
 for XLat = 0, (ModelLatN-1) do begin				; load MapIDLReg, with MissVal or XReg
  for XEight = 0, (LongEightN-1) do begin
    readf, lunRef, Int0, Int1, Int2, Int3, Int4, Int5, Int6, Int7, format="(8(i6))"
    
    MapIDLReg (((XEight*8)+0),XLat) = Int0
    MapIDLReg (((XEight*8)+1),XLat) = Int1
    MapIDLReg (((XEight*8)+2),XLat) = Int2
    MapIDLReg (((XEight*8)+3),XLat) = Int3
    MapIDLReg (((XEight*8)+4),XLat) = Int4
    MapIDLReg (((XEight*8)+5),XLat) = Int5
    MapIDLReg (((XEight*8)+6),XLat) = Int6
    MapIDLReg (((XEight*8)+7),XLat) = Int7
  endfor
 endfor
endif else if ((float(ModelLongN)/10.0) EQ fix(float(ModelLongN)/10.0)) then begin
 LongTenN = round (ModelLongN / 10)
 for XLat = 0, (ModelLatN-1) do begin				; load MapIDLReg, with MissVal or XReg
  for XTen = 0, (LongTenN-1) do begin
    readf, lunRef, Int0, Int1, Int2, Int3, Int4, Int5, Int6, Int7, Int8, Int9, format="(10(i6))"
    
    MapIDLReg (((XTen*8)+0),XLat) = Int0
    MapIDLReg (((XTen*8)+1),XLat) = Int1
    MapIDLReg (((XTen*8)+2),XLat) = Int2
    MapIDLReg (((XTen*8)+3),XLat) = Int3
    MapIDLReg (((XTen*8)+4),XLat) = Int4
    MapIDLReg (((XTen*8)+5),XLat) = Int5
    MapIDLReg (((XTen*8)+6),XLat) = Int6
    MapIDLReg (((XTen*8)+7),XLat) = Int7
    MapIDLReg (((XTen*8)+8),XLat) = Int8
    MapIDLReg (((XTen*8)+9),XLat) = Int9
  endfor
 endfor
endif else if ((float(ModelLongN)/6.0) EQ fix(float(ModelLongN)/6.0)) then begin
 LongTenN = round (ModelLongN / 6)
 for XLat = 0, (ModelLatN-1) do begin				; load MapIDLReg, with MissVal or XReg
  for XTen = 0, (LongTenN-1) do begin
    readf, lunRef, Int0, Int1, Int2, Int3, Int4, Int5, format="(6(i6))"
    
    MapIDLReg (((XTen*6)+0),XLat) = Int0
    MapIDLReg (((XTen*6)+1),XLat) = Int1
    MapIDLReg (((XTen*6)+2),XLat) = Int2
    MapIDLReg (((XTen*6)+3),XLat) = Int3
    MapIDLReg (((XTen*6)+4),XLat) = Int4
    MapIDLReg (((XTen*6)+5),XLat) = Int5
  endfor
 endfor
endif else begin
 print, "  > ##### ERROR: LoadRef: LongN is not divisible by 8 or 10 #####"
endelse

free_lun, lunRef

for XLong = 0, (ModelLongN-1) do begin				; adjust MapIDLReg from 1...RegN
  for XLat = 0, (ModelLatN-1) do begin				;    to 0...(RegN-1)
    if (MapIDLReg (XLong,XLat) NE MissVal) then MapIDLReg (XLong,XLat) = MapIDLReg (XLong,XLat) - 1L
  endfor
endfor

;LogFilePath = '/cru/u2/f709762/data/scratch/log-idl.dat'
;openw, lunLog, LogFilePath, /get_lun			; ###############

TotValid = 0L
ThisReg  = 0L
for XLong = 0, (ModelLongN-1) do begin				; construct RegSize 
  for XLat = 0, (ModelLatN-1) do begin
    ThisReg = long (MapIDLReg (XLong,XLat))
    
;    printf, lunLog, XLong, XLat, ThisReg		; ################
    
    if (ThisReg NE MissVal) then begin
    	RegSize (ThisReg) = RegSize (ThisReg) + 1
    	TotValid          = TotValid          + 1
    endif
  endfor
endfor

;free_lun, lunLog					; ###############

for XReg = 0L, (RegN-1L) do begin				; construct RegColour    
    RegColour (XReg) = 2 + (16 * (XReg MOD 16L))
endfor

print, "  > Regions loaded and total non-missing boxes: ", RegN, TotValid, format='(a48,i8,i8)'

end
