!
!
! Subroutines: Output_Grid, AutoSave
! 
!
subroutine Output_Grid
!
! output grid infor. in the format of .FEG
! input: Linklists of node, element, fault
!
use dflib
use global

implicit none
character(len=512) msg0, msg1
integer(4) iret, ierr
integer(4) i, no, n1, n2, n3, n4
real lon, lat, elev, q, hc,hm, dips1,dips2, offst, density_anomaly, cooling_curvature
real tempvec(3)

msg0 = ' 'C
msg1 = ' 'C
open(unit = 200, file = feg_sav, status = 'new', action='write', iostat = ierr)
if (ierr /= 0) then
  msg0 = 'Specified file already exists. Do you want to overwrite it?'C
  msg1 = 'Error opening file'C
  iret = messageboxqq(msg0, msg1, MB$ICONEXCLAMATION.OR.MB$YESNOCANCEL)
  if (iret == MB$IDYES) then
     open(unit = 200, file = feg_sav, status = 'old', iostat = ierr)
	 rewind(200)
  else
     return
  endif
endif
! write out into file
if (len_trim(title) == 0) then
   title= ' Specify whatever you like here (no more than 80 characters) '
endif
write(200, '(A80)') title
! node
brief =.false.
! consistent with old Orbweave.exe
nrealn = numnod
! if (nrealn == 0) nrealn = numnod
if (n1000 == 0) n1000 = 1000000
!
write(200, 2) numnod, nrealn, nfaken, n1000, brief
2  format(4I8,L8,' (NUMNOD, NREALN, NFAKEN, N1000, BRIEF)')
outputnode: do i = 1, numnod
    tempvec(1) = nodeABG(1,i)
	tempvec(2) = nodeABG(2,i)
	tempvec(3) = nodeABG(3,i)
	call ABG2lonlat(tempvec, lon,lat)
	lon = lon * 57.2957795
	lat = lat * 57.2957795
	elev = eqcm(1,i)
	q = eqcm(2,i)
	hc = eqcm(3,i)
	hm = eqcm(4,i)
    density_anomaly = eqcm(5,i)
    cooling_curvature = eqcm(6,i)
    IF (OrbData5) THEN ! OrbData5 format: 6 real variables per node
	    write(200, 3) i, lon, lat, elev, q, hc, hm, density_anomaly, cooling_curvature
    ELSE ! old OrbData format; only 4 real variables per node
	    write(200, 3) i, lon, lat, elev, q, hc, hm
    END IF
end do outputnode
3  format(I6,0P,2F11.5,1P,6E10.2)
! element
write(200,4) numel
4  format(I10,' (NUMEL = NUMBER OF TRIANGULAR CONTINUUM ELEMENTS)')
outputele: do i = 1, numel
   n1 = nodes(1,i)
   n2 = nodes(2,i)
   n3 = nodes(3,i)
   write(200, '(I6, 3I6)') i, n1, n2, n3
end do outputele
! fault
write(200,5) nfl
5 format (I10,' (NFL =  NUMBER OF CURVILINEAR FAULT ELEMENTS)')
outputfault: do i = 1, nfl
   n1 = nodef(1,i)
   n2 = nodef(2,i)
   n3 = nodef(3,i)
   n4 = nodef(4,i)
   dips1 = fdip(1,i)
   dips2 = fdip(2,i)
   offst = offset(i)
   !
   ! set dip from [1 179] to [-89 89]
   ! but within orbwin, dips are in [1 179], when reading grid, dips in [-89 89]
   ! are changed back to [1 179], see ReadGrid.f90 for details
   !
   if (dips1 > 90.0) dips1 = dips1 - 180.0
   if (dips2 > 90.0) dips2 = dips2 - 180.0
   write(200, '(I6,4I6,2F6.1,1P,E10.2)') i, n1, n2, n3, n4, dips1, dips2, offst
end do outputfault
close(unit = 200)
! success output
msg0 = 'Successfully saved grid as ' //feg_sav//' 'C
msg1 = 'Successfully saved'C
iret = messageboxqq(msg0, msg1, MB$ICONEXCLAMATION.OR.MB$OK)
return
end subroutine Output_Grid
!
!
!
!
subroutine AutoSave(filesave)
use dflib
use global
implicit none
character*(*) filesave
integer(4) iret, ierr
integer(4) i, no, n1, n2, n3, n4
real lon, lat, elev, q, hc,hm, dips1,dips2, offst, density_anomaly, cooling_curvature
real tempvec(3)

open(unit = 400, file = TRIM(filesave), status = 'unknown', action='write', iostat = ierr)
! write out into file
write(400, "(A)") TRIM(filesave)
! node
brief =.false.
! consistent with old Orbweave.exe
nrealn = numnod
! if (nrealn == 0) nrealn = numnod
if (n1000 == 0) n1000 = 100000
!
write(400, 2) numnod, nrealn, nfaken, n1000, brief
2  format(4I8,L8,' (NUMNOD, NREALN, NFAKEN, N1000, BRIEF)')
outputnode: do i = 1, numnod
    tempvec(1) = nodeABG(1,i)
	tempvec(2) = nodeABG(2,i)
	tempvec(3) = nodeABG(3,i)
	call ABG2lonlat(tempvec, lon,lat)
	lon = lon * 57.2957795
	lat = lat * 57.2957795
	elev = eqcm(1,i)
	q = eqcm(2,i)
	hc = eqcm(3,i)
	hm = eqcm(4,i)
    density_anomaly = eqcm(5,i)
    cooling_curvature = eqcm(6,i)
    IF (OrbData5) THEN ! OrbData5 format: 6 real values per node
	    write(400, 3) i, lon, lat, elev, q, hc, hm, density_anomaly, cooling_curvature
    ELSE ! old OrbData format: only 4 real values per node
	    write(400, 3) i, lon, lat, elev, q, hc, hm
    END IF
end do outputnode
3  format(I6,0P,2F11.5,1P,6E10.2)
! element
write(400,4) numel
4  format(I10,' (NUMEL = NUMBER OF TRIANGULAR CONTINUUM ELEMENTS)')
outputele: do i = 1, numel
   n1 = nodes(1,i)
   n2 = nodes(2,i)
   n3 = nodes(3,i)
   write(400, '(I6, 3I6)') i, n1, n2, n3
end do outputele
! fault
write(400,5) nfl
5 format (I10,' (NFL =  NUMBER OF CURVILINEAR FAULT ELEMENTS)')
outputfault: do i = 1, nfl
   n1 = nodef(1,i)
   n2 = nodef(2,i)
   n3 = nodef(3,i)
   n4 = nodef(4,i)
   dips1 = fdip(1,i)
   dips2 = fdip(2,i)
   if (dips1 > 90.0) dips1 = dips1 - 180.0
   if (dips2 > 90.0) dips2 = dips2 - 180.0
   offst = offset(i)
   write(400, '(I6,4I6,2F6.1,1P,E10.2)') i, n1, n2, n3, n4, dips1, dips2, offst
end do outputfault
close(unit = 400)
! success output
msg = ' Successfully saved backup grid in ' // TRIM(filesave)
call setstatusbar(1, TRIM(msg) // CHAR(0))
return
end subroutine Autosave