c--------( GTRF : translate RF-residual between models )-------- parameter (LX=4001, LY=4001) c dimension f(LX,LY) character*80 wdr, flog, fnam1, fnam3 character area*8, narea*8, anul*10, area1*8, salt*8 character srf0*9, srf1*9, buf*80, out*72, sdtm*20 c call premsg('GTRF : translate RF-residual between models') call opnpin() ldr = lwkdir(50,wdr) + 1 flog = wdr fnam1 = wdr fnam3 = wdr call gparma('Enter LOG filename ==> ', 81-ldr, flog(ldr:80)) call gparma('Enter Input source data filename ==> ', * 81-ldr, fnam1(ldr:80)) open(1,file=fnam1,status='old') 1 read(1,'(a)') buf if(buf(1:1).eq.'#') goto 1 read(buf,'(a8,i4,4x,4i8)') area, ncc, iorg,korg, ispa,ispb nc = ncc if(ncc.ge.200) nc = ncc - 200 if(nc.ge.1.and.nc.le.62) then iorg = 0 ispa = 0 ispb = 0 korg = (nc-30)*360 - 180 if(nc.gt.60) korg = 0 if(nc.eq.61) iorg = +5400 if(nc.eq.62) iorg = -5400 endif read(1,'(a)') buf read(buf,*) ixs, iys, mszx, mszy, mx, my, vnul, alt write(salt,'(1x,f7.0)') alt if (alt.eq.0.) salt = ' var.' if (alt.lt.0.) salt = ' undef' xs = float(ixs)/1000. ys = float(iys)/1000. if(mx.gt.LX.or.my.gt.LY) call abendm('too large array size') write(anul,'(f9.2)') vnul if(anul(9:9).ne.'0'.or.anul(1:1).ne.' ') * call abendm('unable to handle vnul value') write(out,'(a,a8,a,i4,4x,4i8)') * '--- Area : ', area, ' Coord.', ncc, iorg,korg, ispa,ispb call premsg(out) write(out,'(a,2a10,2x,2a8,2a6,a7)') * '--- (mesh) ', 'xs', 'ys', 'msz-x', 'msz-y', 'mx', 'my', * 'alt' call premsg(out) write(out,'(a,2f10.2,2i8,2i6,a8)') * '--- source ', xs, ys, mszx, mszy, mx, my, salt call premsg(out) read(1,*) ((f(i,j),i=1,mx),j=1,my) call gparmf('Enter Survey-Year ==> ', year) call premsg('Select OLD RF model Spec. or by Generation') call gparmi(' 1:IGRF(tentative) 2:DGRF 3:PGFR or 0:byGen. ==> ', * mrf0) if (mrf0.lt.0.or.mrf0.gt.3) call abendm('invalid value') srf0 = 'IGRF ' if (mrf0.eq.0) then call gparmi(' Generation Number ==> ', mdl0) if (mdl0.lt.1.or.mdl0.gt.99) call abendm('invalid value') srf0(5:5) = '-' write(srf0(6:7),'(i2)') mdl0 else call gparmi(' RF base year ==> ', mdl0) if (mdl0.lt.1900) call abendm('invalid value') write(srf0(6:9),'(i4)') mdl0 if (mrf0.eq.2) srf0(1:1) = 'D' if (mrf0.eq.3) srf0(1:1) = 'P' endif call premsg('Select NEW RF model Spec. or by Generation') call gparmi(' 1:IGRF(tentative) 2:DGRF 3:PGFR or 0:byGen. ==> ', * mrf1) if (mrf1.lt.0.or.mrf1.gt.3) call abendm('invalid value') srf1 = 'IGRF ' if (mrf1.eq.0) then call gparmi(' Generation Number ==> ', mdl1) if (mdl1.lt.1.or.mdl1.gt.99) call abendm('invalid value') srf1(5:5) = '-' write(srf1(6:7),'(i2)') mdl1 else call gparmi(' RF base year ==> ', mdl1) if (mdl1.lt.1900) call abendm('invalid value') write(srf1(6:9),'(i4)') mdl1 if (mrf1.eq.2) srf0(1:1) = 'D' if (mrf1.eq.3) srf0(1:1) = 'P' endif call gparma('Enter Name-label for new file ==> ', 8, narea) call gparma('Enter output new data filename ==> ', * 81-ldr, fnam3(ldr:80)) call clspin() open(3,file=fnam3,status='new') write(6,'(a)') '----------------------------------------' write(6,'(a)') 'GTRF : translate RF-residual between models' write(6,'(a,a)') ' Source infile : ', fnam1(1:lrtrim(fnam1)) write(6,'(a,a)') ' areaname : ', area write(6,'(a,f8.2)') ' year of survey : ', year write(6,'(a,a8)') ' flying alt. (m) : ', salt write(6,'(a,a)') ' Reference Field : ', srf0 write(6,'(a,a)') ' Output new file : ', fnam3(1:lrtrim(fnam3)) write(6,'(a,a)') ' areaname : ', narea write(6,'(a,a)') ' Reference Field : ', srf1 write(6,1000) narea, ncc, iorg,korg, ispa,ispb, * xs, ys, mszx, mszy, mx, my, salt c write(3,'(a,a)') '# Source: ', fnam1(1:lrtrim(fnam1)) write(3,'(a2,a8,i4,4x,4i8)') * '# ', area, ncc, iorg,korg, ispa,ispb write(3,'(a,f8.2)') '# year of survey : ', year write(3,'(a,a8)') '# flying alt. (m) : ', salt write(3,'(a,a)') '# Reference Field : ', srf0 write(3,'(a)') '# converted into another RF residual' write(3,'(a,a)') '# Reference Field : ', srf1 write(3,'(a8,i4,4x,4i8)') narea, ncc, iorg,korg, ispa,ispb write(3,'(2i12,4i6,1x,f7.1,1x,f7.0)') * ixs,iys, mszx,mszy, mx,my, vnul, alt c call cvinit(ncc, 0.,0.,0.,0.) call cvenik(ys, xs, fi, fk) if (mrf0.eq.0) then call gigrf(mdl0, year) else if (mrf0.eq.1) then call sigrf(float(mdl0)) call tyear(year) else if (mrf0.eq.2) then call sdgrf(float(mdl0)) call tyear(year) else if (mrf0.eq.3) then call spgrf(float(mdl0)) call tyear(year) endif call igrfc(fi/60., fk/60., alt, rdb) call dpcini('... adding old model value: ') do 100 j=1,my y = ys + (j-1)*mszy/1000. do 100 i=1,mx x = xs + (i-1)*mszx/1000. if(f(i,j).ne.vnul) then call cvenik(y, x, fi, fk) call igrfc(fi/60., fk/60., alt, rd) f(i,j) = f(i,j) + (rd-rdb) endif call dpcent((j-1)*mx+i, mx*my) 100 continue if (mrf1.eq.0) then call gigrf(mdl1, year) else if (mrf1.eq.1) then call sigrf(float(mdl1)) call tyear(year) else if (mrf1.eq.2) then call sdgrf(float(mdl1)) call tyear(year) else if (mrf1.eq.3) then call spgrf(float(mdl1)) call tyear(year) endif call dpcini('... subtracting new model : ') do 200 j=1,my y = ys + (j-1)*mszy/1000. do 200 i=1,mx x = xs + (i-1)*mszx/1000. if(f(i,j).ne.vnul) then call cvenik(y, x, fi, fk) call igrfc(fi/60., fk/60., alt, rd) f(i,j) = f(i,j) - (rd-rdb) endif call dpcent((j-1)*mx+i, mx*my) 200 continue do 300 j=1,my write(3,'((f7.1,9(1x,f7.1)))') (f(i,j),i=1,mx) 300 continue kn = 0 c 3 read(1,'(a)',end=8) buf if(buf(1:1).eq.'#') goto 3 read(buf,'(a8,i4,4x,4i8)') area1, ncc1, iorg1,korg1, ispa1,ispb1 nc1 = ncc1 if(ncc1.ge.200) nc1 = ncc1 - 200 if(nc1.ge.1.and.nc1.le.62) then iorg1 = 0 ispa1 = 0 ispb1 = 0 korg1 = (nc1-30)*360 - 180 if(nc1.gt.60) korg1 = 0 if(nc1.eq.61) iorg1 = +5400 if(nc1.eq.62) iorg1 = -5400 endif if(area1.ne.area.or.ncc1.ne.ncc .or. * iorg1.ne.iorg.or.korg1.ne.korg .or. * ispa1.ne.ispa.or.ispb1.ne.ispb) * call abendm('Header1 inconsistent') read(1,*) ixs1, iys1, mszx1, mszy1, mx1, my1, vnul1 if(ixs1.ne.ixs.or.iys1.ne.iys .or. * mszx1.ne.mszx.or.mszy1.ne.mszy .or. * mx1.ne.mx.or.my1.ne.my) * call abendm('Header2 inconsistent') write(anul,'(f10.3)') vnul1 if(anul(9:10).ne.'00'.or.anul(1:1).ne.' ') * call abendm('unable to handle vnul1 value') read(1,*) ((f(i,j),i=1,mx),j=1,my) close(1) write(3,'(a)') '# Var.Obs.Altitude data are copied from Source' write(3,'(a8,i4,4x,4i8)') narea, ncc, iorg,korg, ispa,ispb write(3,'(2i12,4i6,1x,f7.1,1x,f7.0)') * ixs,iys, mszx,mszy, mx,my, vnul1, -1. do 400 j=1,my write(3,'((f7.1,9(1x,f7.1)))') (f(i,j),i=1,mx) 400 continue kn = 1 write(6,'(a)') 'Var.Obs.Altitude data are copied from Source' c 8 close(1) close(3) c if(flog(ldr:ldr).ne.' ') then open(99,file=flog,access='append') write(99,'(//a)') 'GTRF : translate RF-residual between models' write(99,'(a,a)') ' Source infile : ', fnam1(1:lrtrim(fnam1)) write(99,'(a,a)') ' areaname : ', area write(99,'(a,f8.2)') ' year of survey : ', year write(99,'(a,a8)') ' flying alt. (m) : ', salt write(99,'(a,a)') ' Reference Field : ', srf0 write(99,'(a,a)') ' Output new file : ', fnam3(1:lrtrim(fnam3)) write(99,'(a,a)') ' areaname : ', narea write(99,'(a,a)') ' Reference Field : ', srf1 write(99,1000) narea, ncc, iorg,korg, ispa,ispb, * xs, ys, mszx, mszy, mx, my, salt if (kn.eq.1) write(99,'(a)') * 'Var.Obs.Altitude data are copied from Source' call strdtm(sdtm) write(99,'(a,a,a)') '=============== ', sdtm, '===============' close(99) endif stop c 1000 format('Areaname : ', a8, 5x, 'Coord.', i4, 4x, 4i8/ * 3x, ' (mesh) ', ' xs ', ' ys ', * ' msz-x', ' msz-y', ' mx', ' my', ' alt '/ * 3x, 'grid out=in ', 2f10.2, 2i8, 2i6, a8) end