! ! Fortran 95 version of Sib-pair ! ! GUI additions ! Hooks to call JAPI (library calling AWT) ! Java Application Programming Interface 1.0.6 (http://www.japi.de/) ! Main author: Merten Joost ! japi is an open source free software GUI toolkit, which makes it easy ! to develop platform independent applications. Written in JAVA and C, ! provides the JAVA AWT Toolkit to non object oriented Languages incl ! Fortran 77 and Fortran 90 onwards. Development stopped in 2003 ! ! Hooks to call PILIB (library calling GTK+) ! Platform Independent Library for use with Fortran 9x ! PILIB 0.4 (http://www.sourceforge.net/projects/pilib) ! Main author: Martin Hierholzer ! #if JAPI module japi implicit none ! boolean integer,parameter :: j_true = 1 integer,parameter :: j_false = 0 ! alignment integer,parameter :: j_left = 0 integer,parameter :: j_center = 1 integer,parameter :: j_right = 2 integer,parameter :: j_top = 3 integer,parameter :: j_bottom = 4 integer,parameter :: j_topleft = 5 integer,parameter :: j_topright = 6 integer,parameter :: j_bottomleft = 7 integer,parameter :: j_bottomright = 8 ! cursor integer,parameter :: j_default_cursor = 0 integer,parameter :: j_crosshair_cursor = 1 integer,parameter :: j_text_cursor = 2 integer,parameter :: j_wait_cursor = 3 integer,parameter :: j_sw_resize_cursor = 4 integer,parameter :: j_se_resize_cursor = 5 integer,parameter :: j_nw_resize_cursor = 6 integer,parameter :: j_ne_resize_cursor = 7 integer,parameter :: j_n_resize_cursor = 8 integer,parameter :: j_s_resize_cursor = 9 integer,parameter :: j_w_resize_cursor = 10 integer,parameter :: j_e_resize_cursor = 11 integer,parameter :: j_hand_cursor = 12 integer,parameter :: j_move_cursor = 13 ! orientation integer,parameter :: j_horizontal = 0 integer,parameter :: j_vertical = 1 ! fonts integer,parameter :: j_plain = 0 integer,parameter :: j_bold = 1 integer,parameter :: j_italic = 2 integer,parameter :: j_courier = 1 integer,parameter :: j_helvetia = 2 integer,parameter :: j_times = 3 integer,parameter :: j_dialogin = 4 integer,parameter :: j_dialogout = 5 ! colors integer,parameter :: j_black = 0 integer,parameter :: j_white = 1 integer,parameter :: j_red = 2 integer,parameter :: j_green = 3 integer,parameter :: j_blue = 4 integer,parameter :: j_cyan = 5 integer,parameter :: j_magenta = 6 integer,parameter :: j_yellow = 7 integer,parameter :: j_orange = 8 integer,parameter :: j_green_yellow = 9 integer,parameter :: j_green_cyan = 10 integer,parameter :: j_blue_cyan = 11 integer,parameter :: j_blue_magenta = 12 integer,parameter :: j_red_magenta = 13 integer,parameter :: j_dark_gray = 14 integer,parameter :: j_light_gray = 15 integer,parameter :: j_gray = 16 ! borderstyle integer,parameter :: j_none = 0 integer,parameter :: j_linedown = 1 integer,parameter :: j_lineup = 2 integer,parameter :: j_areadown = 3 integer,parameter :: j_areaup = 4 ! mouselistener integer,parameter :: j_moved = 0 integer,parameter :: j_dragged = 1 integer,parameter :: j_pressed = 2 integer,parameter :: j_released = 3 integer,parameter :: j_entererd = 4 integer,parameter :: j_exited = 5 integer,parameter :: j_doubleclick = 6 ! j_moved integer,parameter :: j_resized = 1 integer,parameter :: j_hidden = 2 integer,parameter :: j_shown = 3 ! windowlistener integer,parameter :: j_activated = 0 integer,parameter :: j_deactivated = 1 integer,parameter :: j_opened = 2 integer,parameter :: j_closed = 3 integer,parameter :: j_iconified = 4 integer,parameter :: j_deiconified = 5 integer,parameter :: j_closing = 6 ! imagefileformat integer,parameter :: j_gif = 0 integer,parameter :: j_jpg = 1 integer,parameter :: j_ppm = 2 integer,parameter :: j_bmp = 3 ! ledformat integer,parameter :: j_round = 0 integer,parameter :: j_rect = 1 ! randommax integer,parameter :: j_randmax = 2147483647 ! interface logical,external :: j_start logical,external :: j_connect external :: j_setport external :: j_setdebug integer,external :: j_frame integer,external :: j_button integer,external :: j_graphicbutton integer,external :: j_checkbox integer,external :: j_label integer,external :: j_graphiclabel integer,external :: j_canvas integer,external :: j_panel integer,external :: j_borderpanel integer,external :: j_radiogroup integer,external :: j_radiobutton integer,external :: j_list integer,external :: j_choice integer,external :: j_dialog integer,external :: j_window integer,external :: j_popupmenu integer,external :: j_scrollpane integer,external :: j_hscrollbar integer,external :: j_vscrollbar integer,external :: j_line integer,external :: j_printer integer,external :: j_image external :: j_filedialog external :: j_fileselect integer,external :: j_messagebox integer,external :: j_alertbox integer,external :: j_choicebox2 integer,external :: j_choicebox3 integer,external :: j_progressbar integer,external :: j_led integer,external :: j_sevensegment integer,external :: j_meter external :: j_additem integer,external :: j_textfield integer,external :: j_textarea integer,external :: j_menubar integer,external :: j_menu integer,external :: j_helpmenu integer,external :: j_menuitem integer,external :: j_checkmenuitem external :: j_pack external :: j_print external :: j_playsoundfile external :: j_play integer,external :: j_sound external :: j_setfont external :: j_setfontname external :: j_setfontsize external :: j_setfontstyle external :: j_seperator external :: j_disable external :: j_enable logical,external :: j_getstate integer,external :: j_getrows integer,external :: j_getcolumns integer,external :: j_getselect logical,external :: j_isselect logical,external :: j_isvisible logical,external :: j_isparent logical,external :: j_isresizable external :: j_select external :: j_deselect external :: j_multiplemode external :: j_insert external :: j_remove external :: j_removeitem external :: j_removeall external :: j_setstate external :: j_setrows external :: j_setcolumns external :: j_seticon external :: j_setimage external :: j_setvalue external :: j_setradiogroup external :: j_setunitinc external :: j_setblockinc external :: j_setmin external :: j_setmax external :: j_setdanger external :: j_setslidesize external :: j_setcursor external :: j_setresizable integer,external :: j_getlength integer,external :: j_getvalue integer,external :: j_getdanger integer,external :: j_getscreenheight integer,external :: j_getscreenwidth integer,external :: j_getheight integer,external :: j_getwidth integer,external :: j_getinsets integer,external :: j_getlayoutid integer,external :: j_getinheight integer,external :: j_getinwidth external :: j_gettext external :: j_getitem integer,external :: j_getitemcount external :: j_delete external :: j_replacetext external :: j_appendtext external :: j_inserttext external :: j_settext external :: j_selectall external :: j_selecttext integer,external :: j_getselstart integer,external :: j_getselend external :: j_getseltext integer,external :: j_getcurpos external :: j_setcurpos external :: j_setechochar external :: j_seteditable external :: j_setshortcut external :: j_quit external :: j_kill external :: j_setsize integer,external :: j_getaction integer,external :: j_nextaction external :: j_show external :: j_showpopup external :: j_add external :: j_release external :: j_releaseall external :: j_hide external :: j_dispose external :: j_setpos integer,external :: j_getviewportheight integer,external :: j_getviewportwidth integer,external :: j_getxpos integer,external :: j_getypos external :: j_getpos integer,external :: j_getparentid external :: j_setfocus logical,external :: j_hasfocus integer,external :: j_getstringwidth integer,external :: j_getfontheight integer,external :: j_getfontascent integer,external :: j_keylistener integer,external :: j_getkeycode integer,external :: j_getkeychar integer,external :: j_mouselistener integer,external :: j_getmousex integer,external :: j_getmousey external :: j_getmousepos integer,external :: j_getmousebutton integer,external :: j_focuslistener integer,external :: j_componentlistener integer,external :: j_windowlistener external :: j_setflowlayout external :: j_setborderlayout external :: j_setgridlayout external :: j_setfixlayout external :: j_setnolayout external :: j_setborderpos external :: j_sethgap external :: j_setvgap external :: j_setinsets external :: j_setalign external :: j_setflowfill external :: j_translate external :: j_cliprect external :: j_drawrect external :: j_fillrect external :: j_drawroundrect external :: j_fillroundrect external :: j_drawoval external :: j_filloval external :: j_drawcircle external :: j_fillcircle external :: j_drawarc external :: j_fillarc external :: j_drawline external :: j_drawpolyline external :: j_drawpolygon external :: j_fillpolygon external :: j_drawpixel external :: j_drawstring external :: j_setxor integer,external :: j_getimage external :: j_getimagesource external :: j_drawimagesource integer,external :: j_getscaledimage external :: j_drawimage external :: j_drawscaledimage external :: j_setcolor external :: j_setcolorbg external :: j_setnamedcolor external :: j_setnamedcolorbg integer,external :: j_loadimage logical,external :: j_saveimage external :: j_sync external :: j_beep integer,external :: j_random external :: j_sleep end module japi #endif ! ! time/random number generator seeds ! module rndseed integer :: ix,iy,iz integer :: initix,initiy,initiz end module rndseed ! ! epoch for Julian dates (defaults to 2440588==1970-01-01), ! module julian_epoch double precision :: epoch = 2440588.0d0 end module julian_epoch ! ! interrupt ! module interrupt integer :: irupt contains ! ! keyboard interrupt handler (exits after 6 interrupts, usually ctrl-C) ! subroutine handler() irupt=irupt+1 if (irupt > 5) then write(*,'(a)') 'Multiple interrupts received! Exiting.' stop end if end subroutine handler end module interrupt ! ! Output stream ! module outstream integer :: outstr end module outstream ! ! Input buffer, prompt string etc ! ilevel=0 historical or macro command; =1 keyboard; =2,3 files ! ilevold=last level, if ilevel currently 0 ! module iobuff integer, parameter :: LINSIZ = 20000 character (len=LINSIZ) :: lin, commands character (len=3):: prompt_string = '>> ' integer :: ilevel=1, ilevold=1 end module iobuff ! ! Locus types ! loctyp 1=marker 2=X-marker 3=quant 4=affection 5=mito/Y-marker 11-15=deleted ! module locus_types integer, parameter :: LOC_ANY=0, LOC_DEL=10 integer, parameter :: LOC_CODOM=1, LOC_XLIN=2, LOC_HAP=5, LOC_QUA=3, LOC_AFF=4 integer, parameter :: DEL_CODOM=LOC_DEL+LOC_CODOM, DEL_XLIN=LOC_DEL+LOC_XLIN, & DEL_HAP=LOC_DEL+LOC_HAP, DEL_QUA=LOC_DEL+LOC_QUA, & DEL_AFF=LOC_DEL+LOC_AFF character (len=1), dimension(6), parameter :: typloc = & (/ 'm','x','q','a','h','d' /) character (len=12), dimension(5), parameter :: typlloc(5) = & (/'marker ', 'xmarker ', & 'quantitative', 'affection ', & 'haploid ' /) contains function loccode(ch) integer :: loccode character (len=1), intent(in) :: ch loccode=0 if (ch == 'm') then loccode=LOC_CODOM else if (ch == 'x') then loccode=LOC_XLIN else if (ch == 'q') then loccode=LOC_QUA else if (ch == 'a') then loccode=LOC_AFF else if (ch == 'h') then loccode=LOC_HAP else if (ch == 'd') then loccode=LOC_DEL end if end function loccode end module locus_types ! ! Storage for the locus data ! module locus_data use outstream ! ! Locus structure: ! number of loci, locus name, locus type, and locus position in file ! loctyp 1=marker 2=X-marker 3=quant 4=affection 5=mito/Y-marker 11-15=deleted ! integer :: nloci = 0 character(len=20), dimension(:), allocatable :: loc character(len=40), dimension(:), allocatable :: locnotes integer, dimension(:), allocatable :: loctyp, locpos ! position of locus on sex-averaged linkage map double precision, dimension(:), allocatable :: map integer, dimension(:), allocatable :: locord integer, dimension(:), allocatable :: wloc ! collected statistics, such as P-values double precision, dimension(:), allocatable :: locstat contains ! ! Allocate storage for locus descriptions ! subroutine setup_loci(n) integer :: n integer, parameter :: MISS=-9999 allocate(loc(n), loctyp(n), locpos(n), locnotes(n)) allocate(map(n), locstat(n)) allocate(locord(n), wloc(n)) locstat=MISS end subroutine setup_loci ! ! Deallocate storage for locus descriptions ! subroutine cleanup_loci() deallocate(loc, loctyp, locpos, locnotes) deallocate(map, locstat) deallocate(locord, wloc) end subroutine cleanup_loci ! ! Expand arrays for locus descriptions ! subroutine expand_loci(nextra, plevel) integer, intent(in) :: nextra integer, intent(in) :: plevel integer :: newsiz, oldsiz ! temporary storage character (len=20), dimension(:), allocatable :: loc2 integer, dimension(:), allocatable :: loctyp2 integer, dimension(:), allocatable :: locpos2 character (len=40), dimension(:), allocatable :: locnotes2 double precision, dimension(:), allocatable :: map2 oldsiz=size(loc) newsiz=oldsiz+nextra if (plevel > 0) then write(outstr,*) 'NOTE: expanding maximum number of loci to ', newsiz end if allocate(loc2(oldsiz), loctyp2(oldsiz), locpos2(oldsiz), & locnotes2(oldsiz), map2(oldsiz)) loc2=loc loctyp2=loctyp locpos2=locpos locnotes2=locnotes map2=map call cleanup_loci() call setup_loci(newsiz) loc(1:oldsiz)=loc2(1:oldsiz) loctyp(1:oldsiz)=loctyp2(1:oldsiz) locpos(1:oldsiz)=locpos2(1:oldsiz) locnotes(1:oldsiz)=locnotes2(1:oldsiz) map(1:oldsiz)=map2(1:oldsiz) deallocate(loc2, loctyp2, locpos2, locnotes2, map2) end subroutine expand_loci ! ! insert a locus into the locus list ! other details of the slot are left blank ! subroutine insloc(pos) integer, intent(inout) :: pos integer :: i do i=nloci-1, pos, -1 loc(i+1)=loc(i) loctyp(i+1)=loctyp(i) locpos(i+1)=locpos(i) map(i+1)=map(i) end do end subroutine insloc ! ! Allocate arrays for locus descriptions from MERLIN locus file ! subroutine setupmer(strm) integer, intent(in) :: strm integer :: ioerr, narg, newsiz, oldsiz character (len=100) :: s character (len=40), dimension(2) :: words interface subroutine args(s, narg, arg, typ) character (len=*), intent(in) :: s integer, intent(in) :: typ integer, intent(inout) :: narg character (len=*), dimension(:), intent(out) :: arg end subroutine end interface newsiz=0 do read(strm, '(a)', iostat=ioerr) s if (ioerr /= 0) exit narg=2 call args(s, narg, words, 1) if (words(1)(1:1) /= 'E' .and. words(1)(1:1) /= 'S') then newsiz=newsiz+1 end if end do rewind(strm) oldsiz=size(loc) if (newsiz>oldsiz) then newsiz=5*(1+newsiz/5) write(outstr,*) 'NOTE: expanding maximum number of loci to ', newsiz call expand_loci(newsiz-nloci, 0) end if end subroutine setupmer ! ! initialize locstat (stores test statistic for each locus) ! subroutine setup_stat() integer, parameter :: MISS=-9999 locstat=MISS end subroutine setup_stat end module locus_data ! ! One big pedigree data structure ! Updating size requires copying entire structure ! (hopefully maintaining contiguous storage) ! module idstring_widths integer, parameter :: ped_width = 20 integer, parameter :: id_width = 14 end module idstring_widths module ped_class use idstring_widths type ped_data integer :: nped ! number of pedigrees integer :: nact ! number of active pedigrees integer :: maxsiz ! size of largest pedigree integer :: maxact ! size of largest active pedigree integer :: nobs ! number of records integer :: numloc ! number of columns of locus data integer :: numcol ! number of available columns ! pedigree level data character (len=ped_width), dimension(:), allocatable :: pedigree integer, dimension(:), allocatable :: num integer, dimension(:), allocatable :: nfound integer, dimension(:), allocatable :: actset ! individual level data integer, dimension(:), allocatable :: iped integer, dimension(:), allocatable :: imztwin character (len=id_width), dimension(:), allocatable :: id integer, dimension(:), allocatable :: fa integer, dimension(:), allocatable :: mo integer, dimension(:), allocatable :: sex double precision, dimension(:,:), allocatable :: locus ! useful work arrays -- usually referring to locus being currently analysed logical, dimension(:), allocatable :: untyped end type ped_data contains ! ! allocate pedigree data ! subroutine setup_peds(nped, nobs, numloc, numcol, dataset, astat) integer :: nobs, nped, numloc type (ped_data) :: dataset integer :: astat, iwidth if (allocated(dataset%locus)) then call cleanup_peds(dataset) end if iwidth=max(numcol, numloc) dataset%nped = nped dataset%nact = nped dataset%maxsiz = 0 dataset%maxact = 0 dataset%nobs = nobs dataset%numloc = numloc dataset%numcol = iwidth allocate(dataset%pedigree(nped)) allocate(dataset%num(0:nped)) allocate(dataset%nfound(nped)) allocate(dataset%actset(nped)) dataset%num(0)=0 allocate(dataset%iped(nobs)) allocate(dataset%imztwin(nobs)) allocate(dataset%id(nobs)) allocate(dataset%fa(nobs)) allocate(dataset%mo(nobs)) allocate(dataset%sex(nobs)) allocate(dataset%locus(nobs, iwidth), stat=astat) allocate(dataset%untyped(nobs)) end subroutine setup_peds ! ! copy pedigree data ! subroutine copy_peds(set1, set2) type (ped_data) :: set1, set2 integer :: i set2%nped = set1%nped set2%nact = set1%nact set2%maxsiz = set1%maxsiz set2%maxact = set1%maxact set2%nobs = set1%nobs set2%numloc = set1%numloc set2%numcol = set1%numcol do i=0, set1%nped set2%num(i) = set1%num(i) end do do i=1, set1%nped set2%pedigree(i) = set1%pedigree(i) set2%nfound(i) = set1%nfound(i) set2%actset(i) = set1%actset(i) end do do i=1, set1%nobs set2%iped(i) = set1%iped(i) set2%imztwin(i) = set1%imztwin(i) set2%id(i) = set1%id(i) set2%fa(i) = set1%fa(i) set2%mo(i) = set1%mo(i) set2%sex(i) = set1%sex(i) set2%locus(i, 1:set1%numloc) = set1%locus(i, 1:set1%numloc) end do end subroutine copy_peds ! ! deallocate pedigree structure arrays ! subroutine cleanup_peds(dataset) type (ped_data) :: dataset if (allocated(dataset%locus)) then deallocate(dataset%pedigree) deallocate(dataset%num) deallocate(dataset%nfound) deallocate(dataset%actset) deallocate(dataset%iped) deallocate(dataset%imztwin) deallocate(dataset%id) deallocate(dataset%fa) deallocate(dataset%mo) deallocate(dataset%sex) deallocate(dataset%locus) deallocate(dataset%untyped) end if dataset%nped=0 dataset%nact=0 dataset%maxsiz=0 dataset%maxact=0 dataset%nobs=0 dataset%numloc=0 dataset%numcol=0 end subroutine cleanup_peds ! ! Diagnostics for structure ! subroutine show_ped_allocation(dataset) type (ped_data) :: dataset write(*,*) 'nped =', dataset%nped write(*,*) 'nact =', dataset%nact write(*,*) 'maxsiz=', dataset%maxsiz write(*,*) 'maxact=', dataset%maxact write(*,*) 'nobs =', dataset%nobs write(*,*) 'numloc=', dataset%numloc write(*,*) 'numcol=', dataset%numcol, '=', size(dataset%locus,2) write(*,*) 'Array Allocated?' write(*,*) '-------- ----------' write(*,*) 'pedigree: ', allocated(dataset%pedigree), size(dataset%pedigree) write(*,*) 'num : ', allocated(dataset%num) , size(dataset%num) write(*,*) 'nfound : ', allocated(dataset%nfound) , size(dataset%nfound) write(*,*) 'actset : ', allocated(dataset%actset) , size(dataset%actset) write(*,*) 'iped : ', allocated(dataset%iped) , size(dataset%iped) write(*,*) 'imztwin : ', allocated(dataset%imztwin) , size(dataset%imztwin) write(*,*) 'id : ', allocated(dataset%id) , size(dataset%id) write(*,*) 'fa : ', allocated(dataset%fa) , size(dataset%fa) write(*,*) 'mo : ', allocated(dataset%mo) , size(dataset%mo) write(*,*) 'sex : ', allocated(dataset%sex) , size(dataset%sex) write(*,*) 'locus : ', allocated(dataset%locus) , size(dataset%locus) write(*,*) 'untyped : ', allocated(dataset%untyped) , size(dataset%untyped) end subroutine end module ped_class ! ! Allele frequency data structure ! module alleles_class private public :: allele_data, copyfreq, expand_alleles, cleanup_alleles, genot, calc_gtp_freqs type allele_data integer :: numal ! number of different alleles observed for marker integer :: numgtp ! number of possible genotypes for marker integer :: typed ! number of individuals genotyped at marker integer :: untyped ! number of individuals not genotyped at marker integer :: totall ! number of alleles integer :: topall ! most frequent allele logical :: xlinkd ! sex-linked integer, dimension(:), allocatable :: allele_names ! allele names double precision, dimension(:), allocatable :: allele_freqs ! allele freqs double precision, dimension(:), allocatable :: cum_freqs ! cumulative allele freqs double precision, dimension(:), allocatable :: gtp_freqs ! (log) genotype freqs end type allele_data contains ! ! Copy allele frequency data from one structure to another ! subroutine copyfreq(allele_buffer, allele_buffer2) type (allele_data), intent(inout) :: allele_buffer type (allele_data), intent(inout) :: allele_buffer2 integer :: numal numal=allele_buffer%numal if (allele_buffer2%numal < numal) then deallocate(allele_buffer2%allele_names) deallocate(allele_buffer2%allele_freqs) deallocate(allele_buffer2%cum_freqs) allocate(allele_buffer2%allele_names(numal)) allocate(allele_buffer2%allele_freqs(numal)) allocate(allele_buffer2%cum_freqs(numal)) end if allele_buffer2%allele_names(1:numal)=allele_buffer%allele_names(1:numal) allele_buffer2%allele_freqs(1:numal)=allele_buffer%allele_freqs(1:numal) allele_buffer2%cum_freqs(1:numal)=allele_buffer%cum_freqs(1:numal) allele_buffer2%numal=numal allele_buffer2%numgtp=allele_buffer%numgtp allele_buffer2%topall=allele_buffer%topall allele_buffer2%xlinkd=allele_buffer%xlinkd end subroutine copyfreq ! ! expand size of an allele frequency structure ! subroutine expand_alleles(allele_buffer, nextra) integer, intent(in) :: nextra type (allele_data), intent(inout) :: allele_buffer integer numal type (allele_data) :: allele_buffer2 ! allocate a buffer for old data and copy old data across numal=allele_buffer%numal allocate(allele_buffer2%allele_names(numal)) allocate(allele_buffer2%allele_freqs(numal)) allele_buffer2%allele_names=allele_buffer%allele_names allele_buffer2%allele_freqs=allele_buffer%allele_freqs ! reallocate original structure and bring old data back deallocate(allele_buffer%allele_names) deallocate(allele_buffer%allele_freqs) deallocate(allele_buffer%cum_freqs) allocate(allele_buffer%allele_names(numal+nextra)) allocate(allele_buffer%allele_freqs(numal+nextra)) allocate(allele_buffer%cum_freqs(numal+nextra)) allele_buffer%allele_names(1:numal)=allele_buffer2%allele_names(1:numal) allele_buffer%allele_freqs(1:numal)=allele_buffer2%allele_freqs(1:numal) deallocate(allele_buffer2%allele_names) deallocate(allele_buffer2%allele_freqs) end subroutine expand_alleles ! ! release memory held by an allele frequency structure ! subroutine cleanup_alleles(allele_buffer) type (allele_data), intent(inout) :: allele_buffer if (allele_buffer%numal > 0) then deallocate(allele_buffer%allele_names) deallocate(allele_buffer%allele_freqs) deallocate(allele_buffer%cum_freqs) if (allocated(allele_buffer%gtp_freqs)) then deallocate(allele_buffer%gtp_freqs) end if end if allele_buffer%numal=0 allele_buffer%numgtp=0 allele_buffer%typed=0 allele_buffer%untyped=0 allele_buffer%totall=0 allele_buffer%topall=0 allele_buffer%xlinkd=.false. end subroutine cleanup_alleles ! ! produce genotype frequencies for Metropolis algorithm ! subroutine genot(allele_buffer, gfrq) type (allele_data), intent(in) :: allele_buffer double precision, dimension(allele_buffer%numgtp), intent(out) :: gfrq integer :: i, j, ngtp ngtp=0 do i=1, allele_buffer%numal do j=1, i ngtp=ngtp+1 gfrq(ngtp)=allele_buffer%allele_freqs(i)*allele_buffer%allele_freqs(j) if (i /= j) then gfrq(ngtp)=gfrq(ngtp)+gfrq(ngtp) end if end do end do end subroutine genot ! ! or for sequential imputation subroutine calc_gtp_freqs(allele_buffer) type (allele_data), intent(inout) :: allele_buffer integer :: i, j, ngtp if (allocated(allele_buffer%gtp_freqs)) then deallocate(allele_buffer%gtp_freqs) end if allocate(allele_buffer%gtp_freqs(allele_buffer%numgtp)) ngtp=0 do i=1, allele_buffer%numal do j=1, i ngtp=ngtp+1 allele_buffer%gtp_freqs(ngtp)=allele_buffer%allele_freqs(i)*allele_buffer%allele_freqs(j) if (i /= j) then allele_buffer%gtp_freqs(ngtp)=allele_buffer%gtp_freqs(ngtp)+ & allele_buffer%gtp_freqs(ngtp) end if end do end do end subroutine calc_gtp_freqs end module alleles_class ! ! Hash table for indexing IDs etc ! module idhash_class type hash_table logical :: current = .false. ! is hash known to be up to date integer :: nrec = 0 ! size of table (prime) integer :: primroot = 0 ! constant for "expt hash" probe integer, dimension(:), allocatable :: address end type hash_table contains ! ! Allocate an open addressed hash table ! number of buckets is prime ! probe using exponential hash following Muehlbacher 2004 JUCS 10: 1239-1249 ! subroutine setup_hash(nrec, hashtab, load) integer, intent(in) :: nrec type (hash_table), intent(inout) :: hashtab integer, intent(in) :: load integer, dimension(21), parameter :: tabsizes = (/ & 61, 139, 557, 997, 3023, 4093, 7993, 16381, 49943, 79979, 109943, & 131071, 199999, 399989, 499979, 599999, & 699967, 799999, 899981, 999983, 1099997 /) integer, dimension(21), parameter :: primroots = (/ & 2, 2, 2, 7, 5, 2, 5, 2, 5, 2, 5, & 3, 3, 2, 2, 7, 3, 3, 2, 5, 2 /) integer :: nsiz, tload, w tload=load if (load <= 0) tload=1 if (load >= 100) tload=99 hashtab%current = .false. nsiz=nrec do i=1, 20 if (tabsizes(i) > 100*nsiz/tload) then nsiz=tabsizes(i) w=primroots(i) exit end if end do if (allocated(hashtab%address)) then if (hashtab%nrec < nrec) then deallocate(hashtab%address) allocate(hashtab%address(nsiz)) end if else allocate(hashtab%address(nsiz)) end if hashtab%nrec=nsiz hashtab%primroot=w hashtab%address=0 end subroutine setup_hash ! ! release memory held by a hash table ! subroutine cleanup_hash(hashtab) type (hash_table), intent(inout) :: hashtab if (allocated(hashtab%address)) then deallocate(hashtab%address) end if hashtab%current=.false. hashtab%nrec=0 hashtab%primroot=0 end subroutine cleanup_hash end module idhash_class ! ! List of reserved words ! Token name, left binding power, right binding power, operation ! Pos Name LBP RBP Op (1=unary postfix; 2=binary, infix; 3=if; ! --- ------- --- --- -- 10=zero-arg functions eg rand) ! 0 null 0 0 0 ! 1 ( 200 0 0 ! 2 ) 0 5 0 ! 3 if 0 45 3 ! 4 then 5 25 0 ! 5 else 5 25 0 ! 6 * 120 121 2 ! 7 / 120 121 2 ! 8 + 100 101 2 ! 9 - 100 101 2 ! 10 ^ 139 138 2 ! 11 = 0 0 2 ! 12 not 70 70 1 ! 13 and 65 66 2 ! 14 or 60 61 2 ! 15 < 80 80 2 ! 16 > 80 80 2 ! 17 ge 80 80 2 ! 18 le 80 80 2 ! 19 ne 80 80 2 ! 20 eq 80 80 2 ! 21 neg 138 138 1 ! 22 pos 138 138 1 ! 23 abs 140 140 1 ! 24 sqrt 140 140 1 ! 25 log 140 140 1 ! 26 exp 140 140 1 ! 27 sin 140 140 1 ! 28 cos 140 140 1 ! 29 tan 140 140 1 ! 30 asin 140 140 1 ! 31 acos 140 140 1 ! 32 atan 140 140 1 ! 33 inht 140 140 1 ! 34 int 140 140 1 ! 35 round 140 140 1 ! 36 istyp 140 140 1 ! 37 untyp 140 140 1 ! 38 ishet 140 140 1 ! 39 ishom 140 140 1 ! 40 alla 140 140 1 ! 41 allb 140 140 1 ! 42 rand 0 0 10 ! 43 rnorm 0 0 10 ! 44 pi 0 0 10 ! 45 y 0 0 10 ! 46 n 0 0 10 ! 47 x 0 0 10 ! 48 NUM 0 0 10 ! 49 julian 140 140 1 ! 50 greg 140 140 1 ! 51 log10 140 140 1 ! 52 begin 201 0 0 ! 53 end 0 4 0 ! 54 ; 0 0 0 ! 55 eps 0 0 10 ! 56 pnorm 140 140 1 ! 57 qnorm 140 140 1 ! 58 fact 140 140 1 ! module parser_data integer, parameter :: TOKNUM=58 integer, parameter :: TOK_NULL=0, TOK_LBRACKET=1, TOK_RBRACKET=2, & TOK_IF=3, TOK_THEN=4, TOK_ELSE=5, & TOK_MULT=6, TOK_DIVIDE=7, TOK_ADD=8, & TOK_SUBTRACT=9, TOK_POW=10, TOK_EQUAL=11, & TOK_NOT=12, TOK_AND=13, TOK_OR=14, & TOK_LT=15, TOK_GT=16, TOK_GE=17, & TOK_LE=18, TOK_NE=19, TOK_EQ=20, & TOK_NEG=21, TOK_POS=22, TOK_ABS=23, & TOK_SQRT=24, TOK_LOG=25, TOK_EXP=26, & TOK_SIN=27, TOK_COS=28, TOK_TAN=29, & TOK_ASIN=30, TOK_ACOS=31, TOK_ATAN=32, & TOK_INHT=33, TOK_INT=34, TOK_ROUND=35, & TOK_ISTYP=36, TOK_UNTYP=37, TOK_ISHET=38, & TOK_ISHOM=39, TOK_ALLA=40, TOK_ALLB=41, & TOK_RAND=42, TOK_RNORM=43, TOK_PI=44, & TOK_Y=45, TOK_N=46, TOK_X=47, & TOK_NUM=48, TOK_JULIAN=49, TOK_GREG=50, & TOK_LOG10=51, TOK_BEGIN=52, TOK_END=53, & TOK_COLON=54, TOK_EPS=55, TOK_PNORM=56, & TOK_QNORM=57, TOK_FACT=58 character (len=6), dimension(TOKNUM) :: token = & (/'( ',') ','if ','then ','else ','* ', '/ ', & '+ ','- ','^ ','= ','not ','and ', 'or ', & '< ','> ', & 'ge ','le ','ne ','eq ','neg ','pos ', 'abs ', & 'sqrt ','log ','exp ', & 'sin ','cos ','tan ','asin ','acos ','atan ', 'inht ', & 'int ','round ', & 'istyp ','untyp ','ishet ','ishom ','alla ','allb ', & 'rand ','rnorm ','pi ','y ','n ','x ', 'NUM ', & 'julian','greg ','log10 ','begin ','end ',': ', 'eps ', & 'pnorm ','qnorm ','fact ' /) integer, dimension(0:TOKNUM) :: lbp = & (/0, 200,0,0,5,5,120,120, 100,100,139,0,70,65,60,80,80, & 80,80,80,80,138,138,140,140,140,140, & 140,140,140,140,140,140,140,140,140, & 140,140,140,140,140,140,0,0,0,0,0,0,0, & 140,140,140, 201,0,0,0, 140,140,140 /) integer, dimension(0:TOKNUM) :: rbp = & (/0, 0,5,45,25,25,121,121, 101,101,138,0,70,66,61,80,80, & 80,80,80,80,138,138,140,140,140,140, & 140,140,140,140,140,140,140,140,140, & 140,140,140,140,140,140,0,0,0,0,0,0,0, & 140,140,140, 0,4,0,0, 140,140,140 /) integer, dimension(0:TOKNUM) :: op = & (/0, 0,0,3,0,0,2,2, 2,2,2,2,1,2,2,2,2,2,2,2,2,1,1,1,1,1,1, & 1,1,1,1,1,1,1,1,1,1,1,1,1,1,1, 10,10,10,10,10,10,10, 1,1,1, & 0,0,0,10, 1,1,1 /) ! ! Evaluate type of each term in expression word(farg...larg) and actn ! actn=0 error =1 purely arithmetic =2 legal ! ! Types are: wtyp expr ! ----------------- -------- ! tokens 0 + 0...TOKNUM - ! env 100000 + 1...ENVNUM (value) ! constant 200000 value ! trait data 200000 + 1...NLOCI (value) ! constant genotype 300000 value, value ! genotype data 300000 + 1...NLOCI (value, value) ! MISS 400000 MISS ! missing trait 400000 + 1...NLOCI MISS ! MISS genotype 500000 MISS/MISS ! missing trait 500000 + 1...NLOCI MISS/MISS ! NUM 600000 - ! integer, parameter :: addenv=100000 integer, parameter :: addtra=200000 integer, parameter :: addgen=300000 integer, parameter :: addmtr=400000 integer, parameter :: addmge=500000 integer, parameter :: addnum=600000 ! ! Environmental (automatic) variables for evaluator integer, parameter :: ENVNUM = 13 character (len=6), dimension(ENVNUM) :: env = & (/'female','male ','isfou ','isnon ', & 'num ','nfoun ', & 'anymis', 'anytyp','alltyp','numtyp','famnum','index ','commar' /) end module parser_data ! ! Miniscm ! module scheme_lang use interrupt use outstream use iobuff private public :: cleanup_mem, get_string, get_var, isafun, list_var, init_scheme, repl_scheme integer, parameter :: T_FREE=0, T_STRING=1, T_NUMBER=2, T_SYMBOL=4, T_SYNTAX=8, & T_PROC=16, T_PAIR=32, T_CLOSURE=64, T_CONTINUATION=128, & T_MACRO=256, T_PROMISE=512, T_PORT=2048, T_ATOM=16384, & T_CLRATOM=49151, T_MARK=32768, T_UNMARK=32767, & T_MOVED=-1 integer, parameter :: M_LPAREN=0, M_RPAREN=1, M_DOT=2, M_ATOM=3, M_QUOTE=4, & M_COMMENT=5, M_DQUOTE=6, M_BQUOTE=7, M_COMMA=8, M_ATMARK=9, & M_SHARP=10 integer, parameter :: OP_LOAD=0, OP_T0LVL=1, OP_T1LVL=2, OP_READ=3, OP_VALUEPRINT=4, & OP_EVAL=5, OP_E0ARGS=6, OP_E1ARGS=7, OP_APPLY=8, OP_DOMACRO=9 integer, parameter :: OP_LAMBDA=10, OP_QUOTE=11, OP_DEF0=12, OP_DEF1=13, OP_BEGIN=14, & OP_IF0=15, OP_IF1=16, OP_SET0=17, OP_SET1=18, OP_LET0=19, OP_LET1=20, & OP_LET2=21, OP_LET0AST=22, OP_LET1AST=23, OP_LET2AST=24, OP_LET0REC=25, & OP_LET1REC=26, OP_LET2REC=27, OP_COND0=28, OP_COND1=29, OP_DELAY=30, & OP_AND0=31, OP_AND1=32, OP_OR0=33, OP_OR1=34, OP_C0STREAM=35, OP_C1STREAM=36, & OP_0MACRO=37, OP_1MACRO=38, OP_CASE0=39, OP_CASE1=40, OP_CASE2=41 integer, parameter :: OP_PEVAL=42, OP_PAPPLY=43, OP_CONTINUATION=44, OP_ADD=45, & OP_SUB=46, OP_MUL=47, OP_DIV=48, OP_INTDIV=49, OP_REM=50, OP_MOD=51, & OP_CAR=52, OP_CDR=53, OP_CONS=54, OP_SETCAR=55, & OP_SETCDR=56, OP_NOT=57, OP_BOOL=58, OP_ISINT=59, OP_NULL=60, OP_ZEROP=61, OP_POSP=62, & OP_NEGP=63, OP_NUMEQ=64, OP_LESS=65, OP_GRE=66, OP_LEQ=67, OP_GEQ=68, OP_SYMBOL=69, & OP_NUMBER=70, OP_STRING=71, OP_PROC=72, OP_PAIR=73, OP_EQ=74, OP_EQV=75, & OP_FORCE=76, OP_WRITE=77, OP_DISPLAY=78, OP_NEWLINE=79, OP_ERR0=80, & OP_ERR1=81, OP_REVERSE=82, OP_APPEND=83, OP_PUT=84, OP_GET=85, OP_QUIT=86, & OP_GC=87, OP_GCVERB=88, OP_NEWSEGMENT=89 integer, parameter :: OP_RDSEXPR=90, OP_RDLIST=91, OP_RDDOT=92, OP_RDQUOTE=93, OP_RDQQUOTE=94, & OP_RDUNQUOTE=95, OP_RDUQTSP=96 integer, parameter :: OP_P0LIST=97, OP_P1LIST=98, OP_LIST_LENGTH=99, OP_ASSQ=100, OP_PRINT_WIDTH=101, & OP_P0_WIDTH=102, OP_P1_WIDTH=103, OP_GET_CLOSURE=104, OP_CLOSUREP=105, & OP_MACROP=106 integer, parameter :: OP_EXP=107, OP_LOG=108, OP_SIN=109, OP_COS=110, & OP_TAN=111, OP_ASIN=112, OP_ACOS=113, OP_ATAN=114, & OP_SQRT=115, OP_TRUNCATE=116, OP_ROUND=117, & OP_ABS=118, OP_EXPT=119 integer, parameter :: OP_MIN=120, OP_MAX=121, OP_RANDOM=122 integer, parameter :: OP_MKSTRING=123, OP_STRLEN=124, OP_STRSET=125, & OP_SUBSTR=126, OP_STRAPPEND=127, OP_STRSPLIT=128, & OP_STREQ=129, OP_STRLT=130, OP_STRGT=131, & OP_STRLE=132, OP_STRGE=133, OP_STRFIND=134, & OP_STR2NUM=135, OP_NUM2STR=136, OP_SYM2STR=137, & OP_STR2SYM=138 integer, parameter :: OP_SYSTEM=139, OP_MKPORT=140, OP_CLPORT=141, & OP_RDLINE=142, OP_RUNCMD=143, OP_LSLOCI=144, OP_NLOCI=145, & OP_LOCNAM=146, OP_LOCTYP=147, OP_LOCORD=148, & OP_LOCNOT=149, OP_LOCRANK=150, OP_APROPOS=151, OP_HELP=152 integer :: cell_segment = 500 integer :: currentline=0, eol=0 type mcell integer :: iflag integer :: ivalue integer :: slength character, dimension(:), allocatable :: svalue integer :: keynum integer :: car integer :: cdr end type mcell ! Special addresses integer :: nil=1 ! special cell representing empty cell integer :: t=2 ! special cell representing #t integer :: f=3 ! special cell representing #f integer :: global_env=1 ! pointer to global environment integer :: lambda=1 ! pointer to syntax lambda integer :: quote=1 integer :: qquote=1 integer :: unquote=1 integer :: unquotesp=1 ! registers integer :: scm_args=1 integer :: code=1 integer :: dump=1 integer :: envir=1 ! pointer to symbol table integer :: oblist=1 ! ! Evaluator globals ! integer :: oper=1 ! current operation integer :: tok=1 ! current token integer :: value=1 ! value of current expression integer :: print_flag=1 ! print expression ! ! Memory is memsiz array of mcells: mem(memsiz) ! ! memsiz = current maximum allocatable cells ! nextfree= address of next free cell ! fcell = number of free cells ! integer :: memsiz type (mcell), dimension(:), allocatable :: mem integer :: nextfree integer :: fcells ! ! input streams ! integer, parameter :: MAXPORT=5 integer, dimension(5) :: portaddress = (/21, 22, 23, 24, 25/) integer :: nports = 0 contains ! Memory management subroutine setup_mem(siz) integer, intent(in) :: siz type (mcell), dimension(:), allocatable :: tmpmem integer :: j, slen if (.not.allocated(mem)) then memsiz=siz allocate(mem(memsiz)) do j=1, siz mem(j)%iflag=T_FREE mem(j)%ivalue=0 mem(j)%slength=0 mem(j)%keynum=-1 mem(j)%car=nil mem(j)%cdr=nil end do mem(nil)%iflag=ior(T_ATOM, T_MARK) mem(t)%iflag=ior(T_ATOM, T_MARK) mem(f)%iflag=ior(T_ATOM, T_MARK) nextfree=4 fcells=memsiz-4 else if (siz > memsiz) then allocate(tmpmem(memsiz)) do j=1, memsiz tmpmem(j)%iflag=mem(j)%iflag tmpmem(j)%ivalue=mem(j)%ivalue slen=mem(j)%slength tmpmem(j)%slength=slen if (slen > 0) then allocate(tmpmem(j)%svalue(slen)) tmpmem(j)%svalue=mem(j)%svalue end if tmpmem(j)%keynum=mem(j)%keynum tmpmem(j)%car=mem(j)%car tmpmem(j)%cdr=mem(j)%cdr end do call cleanup_mem() allocate(mem(siz)) do j=1, memsiz mem(j)%iflag=tmpmem(j)%iflag mem(j)%ivalue=tmpmem(j)%ivalue slen=tmpmem(j)%slength mem(j)%slength=slen if (slen > 0) then allocate(mem(j)%svalue(slen)) mem(j)%svalue=tmpmem(j)%svalue deallocate(tmpmem(j)%svalue) end if mem(j)%keynum=tmpmem(j)%keynum mem(j)%car=tmpmem(j)%car mem(j)%cdr=tmpmem(j)%cdr end do do j=memsiz+1, siz mem(j)%iflag=T_FREE mem(j)%ivalue=0 mem(j)%slength=0 mem(j)%car=nil mem(j)%cdr=nil mem(j)%keynum=-1 end do deallocate(tmpmem) if (nextfree == nil) then nextfree=memsiz+1 fcells=siz-memsiz end if memsiz=siz end if end subroutine setup_mem ! ! Clean up memory arrays ! ! Zero a block of cells, including deallocating strings ! subroutine cleanup_bank(sta, fin) integer, intent(in) :: fin, sta integer :: j do j=sta, fin if (allocated(mem(j)%svalue)) then mem(j)%slength=0 deallocate(mem(j)%svalue) else end if mem(j)%iflag=T_FREE mem(j)%ivalue=0 mem(j)%keynum=-1 mem(j)%car=nil mem(j)%cdr=nil end do end subroutine cleanup_bank ! ! Free all memory ! subroutine cleanup_mem() call cleanup_bank(1, memsiz) deallocate(mem) end subroutine cleanup_mem ! ! Mark-sweep garbage collector ! ! Mark cells to be saved ! subroutine gc_mark(a) integer, intent(in) :: a integer :: i, n, p, q, t t=nil p=a 20 continue call setmark(p) if (isatom(p)) goto 60 q=car(p) if (q /= nil .and. .not.ismark(q)) then call setatom(p) call set_car(p, t) t=p p=q goto 20 end if 50 continue q=cdr(p) if (q /= nil .and. .not.ismark(q)) then call set_cdr(p, t) t=p p=q goto 20 end if 60 continue if (t == nil) return q=t if (isatom(q)) then call clratom(q) t=car(q) call set_car(q, p) p=q goto 50 else t=cdr(q) call set_cdr(q, p) p=q goto 60 end if end subroutine gc_mark ! ! Copy all registers to free memory ! Reset pointers from old addresses to new addresses ! subroutine gc(a, b, plevel) integer :: a, b integer, intent(in) :: plevel integer :: i nextfree=1 call gc_mark(nil) call gc_mark(t) call gc_mark(f) call gc_mark(oblist) call gc_mark(global_env) call gc_mark(envir) call gc_mark(scm_args) call gc_mark(code) call gc_mark(dump) call gc_mark(a) call gc_mark(b) call clrmark(nil) fcells=0 nextfree=nil do i=1, memsiz if (ismark(i)) then call clrmark(i) else call cleanup_bank(i, i) call set_cdr(i, nextfree) nextfree = i fcells=fcells+1 end if end do if (plevel > 0) then write(*,*) 'GC recovered ', fcells, ' of ', memsiz, ' cells' end if end subroutine gc ! ! Get next free cell ! function getcell(a, b) integer :: getcell integer, intent(in) :: a, b integer :: newsiz, x if (nextfree == nil) then call gc(a, b, 0) if (nextfree == nil) then newsiz=memsiz+cell_segment call setup_mem(newsiz) call gc(a, b, 0) end if end if x=nextfree nextfree=cdr(x) fcells=fcells-1 call set_car(x, nil) call set_cdr(x, nil) getcell=x end function getcell ! ! Cell type operations ! ! Setting values subroutine set_type(p, iflag) integer, intent(in) :: p integer, intent(in) :: iflag mem(p)%iflag=iflag end subroutine set_type subroutine set_ivalue(p, ivalue) integer, intent(in) :: p integer, intent(in) :: ivalue mem(p)%ivalue=ivalue end subroutine set_ivalue subroutine set_string(p, str) integer, intent(in) :: p character (len=*), intent(in) :: str integer :: i, slen if (allocated(mem(p)%svalue)) then deallocate(mem(p)%svalue) end if slen=len(str) mem(p)%slength=slen allocate(mem(p)%svalue(slen)) do i=1, slen mem(p)%svalue(i)=str(i:i) end do end subroutine set_string subroutine set_substring(p, sta, fin, str) integer, intent(in) :: p integer, intent(in) :: fin, sta character (len=*), intent(in) :: str integer :: i, pos pos=0 do i=sta, min(mem(p)%slength, fin) pos=pos+1 mem(p)%svalue(i)=str(pos:pos) end do end subroutine set_substring ! subroutine set_car(p, icar) integer, intent(in) :: icar, p mem(p)%car=icar end subroutine set_car subroutine set_cdr(p, icdr) integer, intent(in) :: icdr, p mem(p)%cdr=icdr end subroutine set_cdr subroutine set_caar(p, icaar) integer, intent(in) :: icaar, p mem(mem(p)%car)%car=icaar end subroutine set_caar subroutine set_cdar(p, icdar) integer, intent(in) :: icdar, p mem(mem(p)%car)%cdr=icdar end subroutine set_cdar subroutine set_syntaxnum(p, op) integer :: op, p mem(p)%keynum=op end subroutine set_syntaxnum ! ! Getting values ! function get_ivalue(p) integer :: get_ivalue integer, intent(in) :: p get_ivalue=mem(p)%ivalue end function get_ivalue ! ! Strings ! function get_string(p) integer, intent(in) :: p character (len=mem(p)%slength) :: get_string integer :: i get_string=' ' if (mem(p)%slength > 0) then do i=1, mem(p)%slength get_string(i:i)=mem(p)%svalue(i) end do end if end function get_string ! Substring function get_substr(p, sta, fin) integer, intent(in) :: p, sta, fin character (len=(fin-sta)) :: get_substr integer :: i, j get_substr=' ' if (mem(p)%slength > 0) then j=0 do i=sta, fin-1 j=j+1 get_substr(j:j)=mem(p)%svalue(i+1) end do end if end function get_substr ! String length function get_strlen(p) integer :: get_strlen integer, intent(in) :: p get_strlen=mem(p)%slength end function get_strlen ! function get_listlen(p) integer :: get_listlen integer :: p integer :: l, x l=0 x=p do while (ispair(x)) l=l+1 x=cdr(x) end do if (x/=nil) l=-1 get_listlen=l end function get_listlen ! ! allowing access to Scheme environment from Sib-pair ! accessible variables are atomic ! result inserted into passed string ! subroutine get_var(string, pos, fin, istat) character (len=*), intent(inout) :: string integer, intent(inout) :: pos, fin integer, intent(out) :: istat integer :: eos, newlen, reslen, tmp, varend, varsta, x, y character (len=40) :: cbuff eos=len(string) istat=0 varsta=pos+1 varend=fin if (string(varsta:varsta) == '(') varsta=varsta+1 if (string(varend:varend) == ')') varend=varend-1 if (varsta <= varend) then tmp=oblist do while (tmp /= nil) if (ceqstr(string(varsta:varend), caar(tmp))) exit tmp=cdr(tmp) end do if (tmp == nil) then istat=-1 goto 999 end if tmp=car(tmp) x=envir do while (x /= nil) y=car(x) do while (y /= nil) if (caar(y) == tmp) exit y=cdr(y) end do if (y /= nil) exit x=cdr(x) end do if (x /= nil) then if (isstring(cdar(y))) then istat=1 reslen=get_strlen(cdar(y)) else if (isnumber(cdar(y))) then istat=2 reslen=get_strlen(cdar(y)) write(cbuff, '(i20)') get_ivalue(cdar(y)) reslen=len_trim(adjustl(cbuff)) end if newlen=reslen+len_trim(string)-fin+pos-1 if (eos > newlen) then if (istat==1) then string=string(1:(pos-1)) // & get_string(cdar(y)) // & string((fin+1):eos) else string=string(1:(pos-1)) // & trim(adjustl(cbuff)) // & string((fin+1):eos) end if pos=pos+reslen+2 else istat=-2 end if else istat=-1 end if else istat=-3 end if ! ! If failed substitution, close up variable in input string ! 999 continue if (istat < 0) then string=string(1:(pos-1)) // string((fin+1):eos) pos=pos-1 end if end subroutine get_var ! ! List Sib-pair accessible (atomic or pair) variables ! subroutine list_var(typ) integer, intent(in) :: typ integer :: i, obj, s, x, y character (len=1) :: ch character (len=20) :: string obj=oblist do while (obj /= nil) s=car(obj) x=envir do while (x /= nil) y=car(x) do while (y /= nil) if (caar(y) == s) then if (typ==1) then if (isstring(cdar(y))) then write(outstr, '(3a)') get_string(car(s)), '=', get_string(cdar(y)) else if (isnumber(cdar(y))) then write(string, '(i20)') get_ivalue(cdar(y)) write(outstr, '(3a)') get_string(car(s)), '=', adjustl(trim(string)) end if else if (ceqstr('*sp-fun*', car(cdar(y)))) then write(outstr, '(2a/a)', advance='no') & get_string(car(s)), ':', ' ' s=cdr(cdar(y)) do i=1, get_strlen(s) ch=get_substr(s, i-1, i) if (ch == ';') then write(outstr,'(/a)', advance='no') ' ' else write(outstr,'(a)', advance='no') ch end if end do write(outstr,*) end if exit end if y=cdr(y) end do if (y /= nil) exit x=cdr(x) end do obj=cdr(obj) end do end subroutine list_var ! ! Test if a macro function exists ! function isafun(nam) integer :: isafun character (len=*), intent(in) :: nam integer :: tmp, x, y isafun=0 tmp=oblist do while (tmp /= nil) if (ceqstr(trim(nam), caar(tmp))) exit tmp=cdr(tmp) end do if (tmp == nil) then return end if tmp=car(tmp) x=envir do while (x /= nil) y=car(x) do while (y /= nil) if (caar(y) == tmp) exit y=cdr(y) end do if (y /= nil) exit x=cdr(x) end do if (x /= nil) then if (ceqstr('*sp-fun*', car(cdar(y)))) then isafun=cdr(cdar(y)) end if end if end function isafun ! function procnum(p) integer :: procnum integer :: p procnum=mem(p)%ivalue end function procnum function syntaxnum(p) integer :: syntaxnum integer :: p syntaxnum=mem(p)%keynum end function syntaxnum ! ! Underlying primitives for Scheme ! function typeof(p) integer :: typeof integer, intent(in) :: p typeof=mem(p)%iflag end function typeof function isstring(p) logical :: isstring integer, intent(in) :: p isstring=(iand(typeof(p), T_STRING) /= 0) end function isstring function isnumber(p) logical :: isnumber integer, intent(in) :: p isnumber=(iand(typeof(p), T_NUMBER) /= 0) end function isnumber function ispair(p) logical :: ispair integer, intent(in) :: p ispair=(iand(typeof(p), T_PAIR) /= 0) end function ispair ! ! car, cdr etc ! function car(p) integer :: car integer, intent(in) :: p car = mem(p)%car end function car function cdr(p) integer :: cdr integer, intent(in) :: p cdr = mem(p)%cdr end function cdr ! function caar(p) integer :: caar integer, intent(in) :: p caar = mem(p)%car caar = mem(caar)%car end function caar function cadr(p) integer :: cadr integer, intent(in) :: p cadr = mem(p)%cdr cadr = mem(cadr)%car end function cadr function cdar(p) integer :: cdar integer, intent(in) :: p cdar = mem(p)%car cdar = mem(cdar)%cdr end function cdar function cddr(p) integer :: cddr integer, intent(in) :: p cddr = mem(p)%cdr cddr = mem(cddr)%cdr end function cddr function cadar(p) integer :: cadar integer, intent(in) :: p cadar = mem(p)%car cadar = mem(cadar)%cdr cadar = mem(cadar)%car end function cadar function caddr(p) integer :: caddr integer, intent(in) :: p caddr = mem(p)%cdr caddr = mem(caddr)%cdr caddr = mem(caddr)%car end function caddr function cadaar(p) integer :: cadaar integer, intent(in) :: p cadaar = mem(p)%car cadaar = mem(cadaar)%car cadaar = mem(cadaar)%cdr cadaar = mem(cadaar)%car end function cadaar function cadddr(p) integer :: cadddr integer, intent(in) :: p cadddr = mem(p)%cdr cadddr = mem(cadddr)%cdr cadddr = mem(cadddr)%cdr cadddr = mem(cadddr)%car end function cadddr function cddddr(p) integer :: cddddr integer, intent(in) :: p cddddr = mem(p)%cdr cddddr = mem(cddddr)%cdr cddddr = mem(cddddr)%cdr cddddr = mem(cddddr)%cdr end function cddddr ! function issymbol(p) logical :: issymbol integer, intent(in) :: p issymbol=(iand(typeof(p), T_SYMBOL) /= 0) end function issymbol function issyntax(p) logical :: issyntax integer, intent(in) :: p issyntax=(iand(typeof(p), T_SYNTAX) /= 0) end function issyntax function isproc(p) logical :: isproc integer, intent(in) :: p isproc=(iand(typeof(p), T_PROC) /= 0) end function isproc function isclosure(p) logical :: isclosure integer, intent(in) :: p isclosure=(iand(typeof(p), T_CLOSURE) /= 0) end function isclosure function iscontinuation(p) logical :: iscontinuation integer, intent(in) :: p iscontinuation=(iand(typeof(p), T_CONTINUATION) /= 0) end function iscontinuation function ispromise(p) logical :: ispromise integer, intent(in) :: p ispromise=(iand(typeof(p), T_PROMISE) /= 0) end function ispromise function isport(p) logical :: isport integer, intent(in) :: p isport=(iand(typeof(p), T_PORT) /= 0) end function isport ! true or false value functions function istrue(p) logical :: istrue integer :: p istrue=(p /= nil .and. p /= f) end function istrue function isfalse(p) logical :: isfalse integer :: p isfalse=(p /= nil .and. p == f) end function isfalse ! ! Garbage collection ! function isatom(p) logical :: isatom integer, intent(in) :: p isatom=(iand(typeof(p), T_ATOM) /= 0) end function isatom subroutine setatom(p) integer, intent(in) :: p call set_type(p, ior(typeof(p), T_ATOM)) end subroutine setatom subroutine clratom(p) integer, intent(in) :: p call set_type(p, iand(typeof(p), T_CLRATOM)) end subroutine clratom function ismark(p) logical :: ismark integer, intent(in) :: p ismark=(iand(typeof(p), T_MARK) /= 0) end function ismark subroutine setmark(p) integer, intent(in) :: p call set_type(p, ior(typeof(p), T_MARK)) end subroutine setmark subroutine clrmark(p) integer, intent(in) :: p call set_type(p, iand(typeof(p), T_UNMARK)) end subroutine clrmark ! ! Cons ! function cons(reg1, reg2) integer cons integer :: reg1, reg2 cons=getcell(reg1, reg2) call set_type(cons, T_PAIR) mem(cons)%car = reg1 mem(cons)%cdr = reg2 end function cons ! ! Contents of a Lisp string cell equal to a Fortran string ! function ceqstr(cstr, reg) logical :: ceqstr character (len=*) :: cstr integer, intent(in) :: reg integer :: i ceqstr=.true. if (mem(reg)%slength /= len(cstr)) then ceqstr=.false. return else do i=1, len(cstr) if (cstr(i:i) /= mem(reg)%svalue(i)) then ceqstr=.false. return end if end do end if end function ceqstr ! ! Contents of a Lisp string equal to a Lisp string ! function streq(a, b) logical :: streq integer, intent(in) :: a, b integer :: i, slen streq=.true. slen=mem(a)%slength if (slen /= mem(b)%slength) then streq=.false. else do i=1, slen if (mem(a)%svalue(i) /= mem(b)%svalue(i)) then streq=.false. return end if end do end if end function streq ! ! Declare a number ! function mk_number(num) integer :: mk_number integer, intent(in) :: num integer :: tmp tmp=getcell(nil, nil) call set_type(tmp, ior(T_NUMBER, T_ATOM)) call set_ivalue(tmp, num) mk_number=tmp end function mk_number ! ! Declare a string ! function mk_string(nam) integer :: mk_string character (len=*), intent(in) :: nam integer :: tmp tmp=getcell(nil, nil) call set_type(tmp, ior(T_STRING, T_ATOM)) call set_string(tmp, nam) mk_string=tmp end function mk_string ! ! Declare a symbol ! function mk_symbol(nam) integer :: mk_symbol character (len=*), intent(in) :: nam integer :: tmp tmp=oblist do while (tmp /= nil) if (ceqstr(trim(nam), caar(tmp))) then mk_symbol=car(tmp) return end if tmp=cdr(tmp) end do tmp=cons(mk_string(nam), nil) call set_type(tmp, T_SYMBOL) oblist=cons(tmp, oblist) mk_symbol=tmp end function mk_symbol ! ! make symbol or number atom from string ! function mk_atom(str) integer :: mk_atom character (len=*), intent(in) :: str integer :: i, ich, slen logical :: hasdot, isnum ! functions integer :: ival double precision :: fval slen=len(str) hasdot=.false. isnum=.true. do i=1, slen ich = ichar(str(i:i)) if (ich < 48 .or. ich > 57) then if ((ich == 43 .or. ich == 45) .and. i==1 .and. slen>1) then continue else if (ich == 46 .and. .not.hasdot) then hasdot=.true. else isnum=.false. exit end if end if end do if (isnum .and. .not.hasdot) then mk_atom=mk_number(ival(str)) else mk_atom=mk_symbol(str) end if end function mk_atom ! ! Make a constant ! function mk_const(nam) integer :: mk_const character (len=*) :: nam integer :: x if (nam == 't') then mk_const=t else if (nam == 'f') then mk_const=f else if (nam(1:1) == 'd') then read(nam, '(i20)') x mk_const=mk_number(x) end if end function mk_const ! ! make closure, c is code, e is environment ! function mk_closure(c, e) integer :: mk_closure integer :: c, e integer :: x x=getcell(c, e) call set_type(x, T_CLOSURE) call set_car(x, c) call set_cdr(x, e) mk_closure=x end function mk_closure ! ! make continuation ! function mk_continuation(d) integer :: mk_continuation integer :: d integer :: x x=getcell(nil, d) call set_type(x, T_CONTINUATION) call set_cdr(x, d) mk_continuation=x end function mk_continuation ! ! make a port ! function mk_port(iport, nam) integer :: mk_port integer, intent(in) :: iport character (len=*), intent(in) :: nam integer :: tmp tmp=getcell(nil, nil) call set_type(tmp, ior(T_PORT, T_ATOM)) call set_ivalue(tmp, iport) call set_string(tmp, nam) mk_port=tmp end function mk_port ! ! Test a port - return location in portaddress function test_port(iport) integer :: test_port integer, intent(in) :: iport integer :: i do i=1, nports if (portaddress(i) == iport) then test_port=i return end if end do test_port=0 end function test_port ! ! Reverse list ! function reverse(a) integer :: reverse integer, intent(in) :: a integer :: p, tmp integer :: i p=nil tmp=a i=0 do while (ispair(tmp)) i=i+1 p=cons(car(tmp), p) tmp=cdr(tmp) end do reverse=p end function reverse ! ! Reverse list -- no new cell generated ! function non_alloc_rev(term, list) integer :: non_alloc_rev integer, intent(in) :: term, list integer :: i, p, res, q i=0 p=list res=term do while (p /= nil) i=i+1 q=cdr(p) call set_cdr(p, res) res=p p=q end do non_alloc_rev=res end function non_alloc_rev ! ! append list -- make new cells ! function append(a, b) integer :: append integer, intent(in) :: a, b integer :: p, q, tmp p=b tmp=a if (tmp /= nil) then tmp = reverse(tmp) do while (tmp /= nil) q = cdr(tmp) call set_cdr(tmp, p) p = tmp tmp = q end do end if append=p end function append ! ! equivalence of atoms ! function eqv(a, b) logical :: eqv integer, intent(in) :: a, b eqv=.false. if (isstring(a)) then if (isstring(b)) then eqv=streq(a, b) end if else if (isnumber(a)) then if (isnumber(b)) then eqv=(get_ivalue(a) == get_ivalue(b)) end if else eqv=(a == b) end if end function eqv ! ! get a new character from input file or stdin ! subroutine inchar(ch) character (len=1) :: ch integer :: ios if (eol==0 .or. currentline > eol) then read(*, '(a)', iostat=ios) lin if (ios /= 0) then write(*, '(a)') 'Exiting!' return end if currentline=0 eol=len_trim(lin) end if currentline=currentline+1 if (currentline > eol) then ch=' ' else ch=lin(currentline:currentline) end if end subroutine inchar ! ! clear input buffer ! subroutine clearinput() currentline=eol end subroutine clearinput ! ! back to standard input ! subroutine flushinput() call clearinput() end subroutine flushinput ! ! backstep one character in input buffer ! subroutine backchar() currentline=currentline-1 end subroutine backchar ! ! skip whitespace ! subroutine skipspace() do currentline=currentline+1 if (currentline > eol) exit if (lin(currentline:currentline)/=' ' .or. & lin(currentline:currentline)/='\t') then exit end if end do currentline=currentline-1 end subroutine skipspace ! ! get next token ! function token() integer :: token character (len=1) :: ch integer :: ich do call skipspace() call inchar(ch) if (ch /= ' ' .and. ch /= '\t') exit end do ich = ichar(ch) if (ch == '(') then token=M_LPAREN else if (ch == ')') then token=M_RPAREN else if (ch == '.') then token=M_DOT else if (ich == 39) then token=M_QUOTE else if (ch == ';') then token=M_COMMENT else if (ch == '"') then token=M_DQUOTE else if (ch == '`') then token=M_BQUOTE else if (ch == ',') then call inchar(ch) if (ch == '@') then token=M_ATMARK else call backchar() token=M_COMMA end if else if (ch == '#') then token=M_SHARP else call backchar() token=M_ATOM end if end function token ! ! read characters to delimiter -- hard coded to work on Windows as well ! function scheme_delim(ch) logical :: scheme_delim character(len=1) :: ch integer :: ich ! functions integer :: ichar ich=ichar(ch) scheme_delim=(ich == 9 .or. ich == 10 .or. & ich == 32 .or. ich== 40 .or. ich == 41) end function scheme_delim ! subroutine readstr(res) character (len=*) :: res integer :: i, pos, reslen character (len=1) :: ch reslen=len(res) res=' ' pos=0 rdloop: do call inchar(ch) if (currentline > eol) exit rdloop if (scheme_delim(ch)) exit rdloop pos=pos+1 if (pos <= reslen) res(pos:pos)=ch end do rdloop call backchar() end subroutine readstr ! ! read rest of a quoted string ! subroutine readstrexp(res, reslen) character (len=*) :: res character (len=1) :: ch integer, intent(out) :: reslen integer :: pos reslen=len(res) res=' ' pos=0 do call inchar(ch) #ifndef IFORT if (ch == '\\') then #else if (ch == '\') then #endif call inchar(ch) else if (ch == '"') then exit end if pos=pos+1 if (pos <= reslen) res(pos:pos)=ch end do reslen=min(pos, reslen) end subroutine readstrexp ! ! print an atom ! subroutine printatom(l) integer :: l character (len=1000) :: str if (l == nil) then write(outstr, '(a)', advance='no') '()' else if (l == t) then write(outstr, '(a)', advance='no') '#t' else if (l == f) then write(outstr, '(a)', advance='no') '#f' else if (isnumber(l)) then write(str, '(i20)') get_ivalue(l) write(outstr, '(a)', advance='no') trim(adjustl(str)) else if (isstring(l)) then if (print_flag==1) then write(outstr,'(3a)',advance='no') '"', get_string(l), '"' else write(outstr,'(a)',advance='no') get_string(l) end if else if (issymbol(l)) then write(outstr,'(a)',advance='no') get_string(car(l)) // ' ' else if (isproc(l)) then write(outstr,'(a,i4,a)', advance='no') '#' else if (isclosure(l)) then write(outstr,'(a)',advance='no') '#' else if (iscontinuation(l)) then write(outstr,'(a)',advance='no') '#' else if (isport(l)) then write(outstr,'(a,i4,3a)',advance='no') '#' end if end subroutine printatom ! function ok_abbrev(x) logical :: ok_abbrev integer :: x ok_abbrev=(ispair(x) .and. cdr(x) /= nil) end function ok_abbrev ! subroutine s_save(a, b, c) integer :: a, b, c dump = cons(envir, cons(c, dump)) dump = cons(b, dump) dump = cons(mk_number(a), dump) end subroutine s_save subroutine s_return(a) integer :: a value = a oper = get_ivalue(car(dump)) scm_args = cadr(dump) envir = caddr(dump) code = cadddr(dump) dump = cddddr(dump) end subroutine s_return subroutine s_retbool(tf) logical :: tf if (tf) then call s_return(t) else call s_return(f) end if end subroutine s_retbool ! ! Apply Scheme commands - split into opexe0 to opexe10 ! subroutine opexe0(op, inline) integer :: op integer :: inline integer :: x, y logical :: ios ! top level if (op == OP_T0LVL) then if (inline /= 3) write(outstr,*) dump = nil envir = global_env call s_save(OP_VALUEPRINT, nil, nil); call s_save(OP_T1LVL, nil, nil); if (inline == 1) write(*, '(a)', advance='no') prompt_string oper=OP_READ else if (op == OP_T1LVL) then code=value oper=OP_EVAL ! read else if (op == OP_READ) then tok=token() oper=OP_RDSEXPR ! print evaluation result else if (op == OP_VALUEPRINT) then print_flag=1 scm_args=value call s_save(OP_T0LVL, nil, nil) oper=OP_P0LIST ! main part of evaluation else if (op == OP_EVAL) then ! symbol if (issymbol(code)) then x=envir do while (x /= nil) y=car(x) do while (y /= nil) if (caar(y) == code) exit y=cdr(y) end do if (y /= nil) exit x=cdr(x) end do if (x /= nil) then call s_return(cdar(y)) else call error1('Unbound variable', code) end if else if (ispair(code)) then ! syntax x = car(code) if (issyntax(x)) then code = cdr(code); oper=syntaxnum(x) ! first, eval top element and eval arguments else call s_save(OP_E0ARGS, nil, code) code = car(code) oper=OP_EVAL end if else call s_return(code) end if ! eval arguments else if (op == OP_E0ARGS) then code=cdr(code) oper=OP_E1ARGS else if (op == OP_E1ARGS) then scm_args=cons(value, scm_args) ! continue if (ispair(code)) then call s_save(OP_E1ARGS, scm_args, cdr(code)) code=car(code) scm_args=nil oper=OP_EVAL ! end else scm_args=reverse(scm_args) code=car(scm_args) scm_args=cdr(scm_args) oper=OP_APPLY end if ! apply code to args else if (op == OP_APPLY) then if (isproc(code)) then oper=procnum(code) else if (isclosure(code)) then envir= cons(nil, cdr(code)) x=car(car(code)) y=scm_args do while (ispair(x)) if (y == nil) then call error0('Too few arguments') return else call set_car(envir, cons(cons(car(x), car(y)), car(envir))) end if x=cdr(x) y=cdr(y) end do if (x == nil) then continue else if (issymbol(x)) then call set_car(envir, cons(cons(x,y), car(envir))) else call error0('Syntax error in closure') return end if code=cdar(code) scm_args=nil oper=OP_BEGIN ! continuation else if (iscontinuation(code)) then dump=cdr(code) if (scm_args /= nil) then call s_return(car(scm_args)) else call s_return(nil) end if else call error0('Illegal function') return end if ! lambda else if (op == OP_LAMBDA) then call s_return(mk_closure(code, envir)) ! quote else if (op == OP_QUOTE) then call s_return(car(code)) ! define else if (op == OP_DEF0) then if (ispair(car(code))) then x=caar(code) code=cons(lambda, cons(cdar(code), cdr(code))) else x=car(code) code=cadr(code) end if if (.not.issymbol(x)) then call error0('Variable is not symbol') return end if call s_save(OP_DEF1, nil, x) oper=OP_EVAL ! define else if (op == OP_DEF1) then x=car(envir) do while (x /= nil) if (caar(x) == code) exit x=cdr(x) end do if (x /= nil) then call set_cdar(x, value) else call set_car(envir, cons(cons(code, value), car(envir))) end if call s_return(code) ! set! else if (op == OP_SET0) then call s_save(OP_SET1, nil, car(code)) code=cadr(code) oper=OP_EVAL ! set! else if (op == OP_SET1) then x=envir do while (x /= nil) y=car(x) do while (y /= nil) if (caar(y) == code) exit y=cdr(y) end do if (y /= nil) exit x=cdr(x) end do if (x /= nil) then call set_cdar(y, value) call s_return(value) else call error1('Unbound variable', code) end if ! begin else if (op == OP_BEGIN) then if (.not.ispair(code)) then call s_return(code) end if if (cdr(code) /= nil) then call s_save(OP_BEGIN, nil, cdr(code)) end if code = car(code) oper=OP_EVAL ! if else if (op == OP_IF0) then call s_save(OP_IF1, nil, cdr(code)) code = car(code) oper=OP_EVAL else if (op == OP_IF1) then if (istrue(value)) then code=car(code) else code=cadr(code) end if oper=OP_EVAL ! let else if (op == OP_LET0) then scm_args=nil value=code if (issymbol(car(code))) then code=cadr(code) else code=car(code) end if oper=OP_LET1 ! let (calculate parameters) else if (op == OP_LET1) then scm_args=cons(value, scm_args) if (ispair(code)) then call s_save(OP_LET1, scm_args, cdr(code)) code=cadar(code) scm_args=nil oper=OP_EVAL else scm_args=reverse(scm_args) code=car(scm_args) scm_args=cdr(scm_args) oper=OP_LET2 end if else if (op == OP_LET2) then envir=cons(nil, envir) if (issymbol(car(code))) then x=cadr(code) else x=car(code) end if y=scm_args do while (y /= nil) call set_car(envir, cons(cons(caar(x), car(y)), car(envir))) x=cdr(x) y=cdr(y) end do ! named let if (issymbol(car(code))) then x=cadr(code) scm_args=nil do while (x /= nil) scm_args=cons(caar(x), scm_args) x=cdr(x) end do x=mk_closure(cons(reverse(scm_args), cddr(code)), envir) call set_car(envir, cons(cons(car(code), x), car(envir))) code=cddr(code) scm_args=nil else code=cdr(code) scm_args=nil end if oper=OP_BEGIN ! let* else if (op == OP_LET0AST) then if (car(code) == nil) then envir = cons(nil, envir) code=cdr(code) oper=OP_BEGIN end if call s_save(OP_LET1AST, cdr(code), car(code)) code=cadaar(code) oper=OP_EVAL ! let* (make new frame) else if (op == OP_LET1AST) then envir=cons(nil, envir) oper=OP_LET2AST ! let* (calculate parameters) else if (op == OP_LET2AST) then call set_car(envir, cons(cons(caar(code), value), car(envir))) code = cdr(code) ! continue if (ispair(code)) then call s_save(OP_LET2AST, scm_args, code) code = cadar(code) scm_args = nil oper=OP_EVAL ! end else code = scm_args scm_args = nil oper=OP_BEGIN end if else write(*, '(i3,a)') oper, ' is an illegal operator' end if end subroutine opexe0 subroutine opexe1(op) integer :: op integer :: x, y ! letrec if (op == OP_LET0REC) then envir = cons(nil, envir) scm_args = nil value = code code = car(code) oper=OP_LET1REC ! letrec calculate parameters else if (op == OP_LET1REC) then scm_args = cons(value, scm_args) if (ispair(code)) then ! continue call s_save(OP_LET1REC, scm_args, cdr(code)) code=cadar(code) scm_args=nil oper=OP_EVAL else ! end scm_args = reverse(scm_args) code = car(scm_args) scm_args = cdr(scm_args) oper=OP_LET2REC end if ! letrec else if (op == OP_LET2REC) then x=car(code) y=scm_args do while (y /= nil) call set_car(envir, cons(cons(caar(x), car(y)), car(envir))) x=cdr(x) y=cdr(y) end do code=cdr(code) scm_args=nil oper=OP_BEGIN ! cond else if (op == OP_COND0) then if (.not.ispair(code)) then call error0('Syntax error in cond') return end if call s_save(OP_COND1, nil, code) code=caar(code) oper=OP_EVAL else if (op == OP_COND1) then if (istrue(value)) then code= cdar(code) if (code == nil) then call s_return(value) end if oper=OP_BEGIN else code=cdr(code) if (code == nil) then call s_return(nil) else call s_save(OP_COND1, nil, code) code=caar(code) oper=OP_EVAL end if end if ! delay else if (op == OP_DELAY) then x=mk_closure(cons(nil, code), envir) call set_type(x, ior(T_PROMISE, typeof(x))) call s_return(x) ! and else if (op == OP_AND0) then if (code==nil) then call s_return(t) end if call s_save(OP_AND1, nil, cdr(code)) code=car(code) oper=OP_EVAL ! and else if (op == OP_AND1) then if (isfalse(value)) then call s_return(value) else if (code == nil) then call s_return(value) else call s_save(OP_AND1, nil, cdr(code)) code=car(code) oper=OP_EVAL end if ! or else if (op == OP_OR0) then if (code==nil) then call s_return(f) end if call s_save(OP_OR1, nil, cdr(code)) code=car(code) oper=OP_EVAL ! or else if (op == OP_OR1) then if (istrue(value)) then call s_return(value) else if (code == nil) then call s_return(value) else call s_save(OP_OR1, nil, cdr(code)) code=car(code) oper=OP_EVAL end if ! cons-stream else if (op == OP_C0STREAM) then call s_save(OP_C1STREAM, nil, cdr(code)) code=car(code) oper=OP_EVAL ! cons-stream else if (op == OP_C1STREAM) then scm_args=value x=mk_closure(cons(nil, code), envir) call set_type(x, ior(T_PROMISE, typeof(x))) call s_return(cons(scm_args, x)) ! case else if (op == OP_CASE0) then call s_save(OP_CASE1, nil, cdr(code)) code=car(code) oper=OP_EVAL ! case else if (op == OP_CASE1) then x=code do while (x /= nil) y=caar(x) if (.not.ispair(y)) then exit end if do while (y /= nil) if (eqv(car(y), value)) then exit end if y=cdr(y) end do if (y /= nil) exit x = cdr(x) end do if (x /= nil) then if (ispair(caar(x))) then code = cdar(x) oper=OP_BEGIN ! else else call s_save(OP_CASE2, nil, cdar(x)) code=caar(x) oper=OP_EVAL end if else call s_return(nil) end if ! case else if (op == OP_CASE2) then if (istrue(value)) then oper=OP_BEGIN else call s_return(nil) end if ! apply else if (op == OP_PAPPLY) then code=car(scm_args) scm_args=cadr(scm_args) oper=OP_APPLY ! eval else if (op == OP_PEVAL) then code=car(scm_args) scm_args=nil oper=OP_EVAL ! call-with-current-continuation else if (op == OP_CONTINUATION) then code=car(scm_args) scm_args=cons(mk_continuation(dump), nil) oper=OP_APPLY else write(outstr, '(2a)') oper, 'is an illegal operator!' end if end subroutine opexe1 ! subroutine opexe2(op) integer :: op integer :: d, x, v integer :: i logical :: int_op double precision :: rv int_op=.true. ! + if (op == OP_ADD) then x = scm_args v = 0 do while (x /= nil) if (.not.isnumber(car(x))) then call error1('Argument is not a number: ', car(x)) return else v=v+get_ivalue(car(x)) end if x=cdr(x) end do call s_return(mk_number(v)) ! - else if (op == OP_SUB) then x = cdr(scm_args) if (.not.isnumber(car(scm_args))) then call error1('Argument is not a number: ', car(scm_args)) return else v = get_ivalue(car(scm_args)) end if i = 1 do while (x /= nil) i=i+1 if (.not.isnumber(car(x))) then call error1('Argument is not a number: ', car(x)) return else v=v-get_ivalue(car(x)) end if x=cdr(x) end do if (i==1) then v=-v end if call s_return(mk_number(v)) ! * else if (op == OP_MUL) then x = scm_args v = 1 do while (x /= nil) if (.not.isnumber(car(x))) then call error1('Argument is not a number: ', car(x)) return else v=v*get_ivalue(car(x)) end if x=cdr(x) end do call s_return(mk_number(v)) ! / quotient else if (op == OP_DIV .or. op == OP_INTDIV) then x = cdr(scm_args) if (.not.isnumber(car(scm_args))) then call error1('Argument is not a number: ', car(scm_args)) return else v = get_ivalue(car(scm_args)) end if do while (x /= nil) if (.not.isnumber(car(x))) then call error1('Argument is not a number: ', car(x)) return else d=get_ivalue(car(x)) end if if (d /= 0) then v=v/d else call error0('Divided by zero!') return end if x=cdr(x) end do call s_return(mk_number(v)) ! remainder else if (op == OP_REM) then x = cdr(scm_args) if (.not.isnumber(car(scm_args))) then call error1('Argument is not a number: ', car(scm_args)) return else v = get_ivalue(car(scm_args)) end if do while (x /= nil) if (.not.isnumber(car(x))) then call error1('Argument is not a number: ', car(x)) return else if (get_ivalue(car(x)) /= 0) then v=mod(v,get_ivalue(car(x))) else call error0('Divided by zero!') return end if end if x=cdr(x) end do call s_return(mk_number(v)) ! modulo else if (op == OP_MOD) then if (isnumber(car(scm_args)) .and. isnumber(cadr(scm_args))) then i = get_ivalue(cadr(scm_args)) if (i /= 0) then v = mod(get_ivalue(car(scm_args)), i) if (v*i < 0) then i=abs(i) if (v > 0) then v=v-i else v=v+i end if end if call s_return(mk_number(v)) else call error0('Modulo x 0 not allowed!') end if end if ! car else if (op == OP_CAR) then if (ispair(car(scm_args))) then call s_return(caar(scm_args)) else call error0('Unable to car for a non-cons cell!') end if ! cdr else if (op == OP_CDR) then if (ispair(car(scm_args))) then call s_return(cdar(scm_args)) else call error0('Unable to cdr for a non-cons cell!') end if ! cons else if (op == OP_CONS) then call set_cdr(scm_args, cadr(scm_args)) call s_return(scm_args) ! set-car! else if (op == OP_SETCAR) then if (ispair(car(scm_args))) then call set_caar(scm_args, cadr(scm_args)) call s_return(car(scm_args)) else call error0('Unable to set-car! for a non-cons cell!') end if ! set-cdr! else if (op == OP_SETCDR) then if (ispair(car(scm_args))) then call set_cdar(scm_args, cadr(scm_args)) call s_return(car(scm_args)) else call error0('Unable to set-cdr! for a non-cons cell!') end if else write(outstr, '(2a)') oper, 'is an illegal operator!' end if end subroutine opexe2 ! subroutine opexe3(op) integer :: op integer :: x, v logical :: comp ! not if (op == OP_NOT) then call s_retbool(isfalse(car(scm_args))) ! boolean? else if (op == OP_BOOL) then call s_retbool(car(scm_args) == f .or. car(scm_args) == t) ! integer? else if (op == OP_ISINT) then call s_retbool(isnumber(car(scm_args))) ! null else if (op == OP_NULL) then call s_retbool(car(scm_args) == nil) ! zero? else if (op == OP_ZEROP) then call s_retbool(get_ivalue(car(scm_args)) == 0) ! positive? else if (op == OP_POSP) then call s_retbool(get_ivalue(car(scm_args)) > 0) ! negative? else if (op == OP_NEGP) then call s_retbool(get_ivalue(car(scm_args)) < 0) ! =, <, >, <=, >= else if (op >= OP_NUMEQ .and. op <= OP_GEQ) then x = cdr(scm_args) if (.not.isnumber(car(scm_args))) then call error1('Argument is not a number: ', car(scm_args)) return else v = get_ivalue(car(scm_args)) end if do while (x /= nil) if (.not.isnumber(car(x))) then call error1('Argument is not a number: ', car(x)) return else if (op == OP_NUMEQ) then comp=(v==get_ivalue(car(x))) else if (op == OP_LESS) then comp=(vget_ivalue(car(x))) else if (op == OP_LEQ) then comp=(v<=get_ivalue(car(x))) else if (op == OP_GEQ) then comp=(v>=get_ivalue(car(x))) end if if (.not.comp) exit end if x=cdr(x) end do call s_retbool(comp) ! symbol? else if (op == OP_SYMBOL) then call s_retbool(issymbol(car(scm_args))) ! number? else if (op == OP_NUMBER) then call s_retbool(isnumber(car(scm_args))) ! string? else if (op == OP_STRING) then call s_retbool(isstring(car(scm_args))) ! procedure? else if (op == OP_PROC) then call s_retbool(isproc(car(scm_args)) .or. isclosure(car(scm_args)) .or. & iscontinuation(car(scm_args))) ! pair? else if (op == OP_PAIR) then call s_retbool(ispair(car(scm_args))) ! eq? else if (op == OP_EQ) then call s_retbool(car(scm_args) == cadr(scm_args)) ! eqv? else if (op == OP_EQV) then call s_retbool(eqv(car(scm_args),cadr(scm_args))) else write(outstr, '(2a)') oper, 'is an illegal operator!' end if end subroutine opexe3 ! subroutine opexe4(op, plevel) integer, intent(in) :: op integer, intent(in) :: plevel integer :: x, y ! force if (op == OP_FORCE) then code = car(scm_args) if (ispromise(code)) then scm_args=nil oper=OP_APPLY else call s_return(code) end if ! write else if (op == OP_WRITE) then print_flag=1 scm_args=car(scm_args) oper=OP_P0LIST ! display else if (op == OP_DISPLAY) then print_flag=0 scm_args=car(scm_args) oper=OP_P0LIST ! newline else if (op == OP_NEWLINE) then if (plevel >= 0) write(outstr,*) call s_return(t) ! error else if (op == OP_ERR0) then if (.not.isstring(car(scm_args))) then call error0('Error -- first argument must be a string') end if write(outstr, '(2a)', advance='no') 'Error: ', get_string(car(scm_args)) scm_args=cdr(scm_args) oper=OP_ERR1 ! error else if (op == OP_ERR1) then write(outstr, '(a)', advance='no') ' ' if (scm_args /= nil) then call s_save(OP_ERR1, cdr(scm_args), nil) scm_args=car(scm_args) print_flag=1 oper=OP_P0LIST else write(outstr,*) call flushinput() oper=OP_T0LVL end if ! reverse else if (op == OP_REVERSE) then call s_return(reverse(car(scm_args))) ! append else if (op == OP_APPEND) then call s_return(append(car(scm_args), cadr(scm_args))) ! gc else if (op == OP_GC) then call gc(nil, nil, 1) call s_return(t) ! new segment size else if (op == OP_NEWSEGMENT) then if (.not.isnumber(car(scm_args))) then call error0('new segment -- argument must be a number!') return end if cell_segment=get_ivalue(car(scm_args)) write(outstr, '(a,i4,a)') & 'Allocation of new memory in increments of ', cell_segment, ' cells' call s_return(t) end if end subroutine opexe4 ! subroutine opexe5(op, plevel) integer, intent(in) :: op integer, intent(in) :: plevel integer :: x integer :: strlen character (len=1000) :: str character (len=1) :: ch if (op == OP_RDSEXPR) then if (tok == M_COMMENT) then do call inchar(ch) if (currentline > eol) exit end do tok = token() oper=OP_RDSEXPR else if (tok == M_LPAREN) then tok = token() if (tok == M_RPAREN) then call s_return(nil) else if (tok == M_DOT) then call error0('Syntax error -- illegal dot expression') else call s_save(OP_RDLIST, nil, nil) oper=OP_RDSEXPR end if else if (tok == M_QUOTE) then call s_save(OP_RDQUOTE, nil, nil) tok=token() oper=OP_RDSEXPR else if (tok == M_BQUOTE) then call s_save(OP_RDQQUOTE, nil, nil) tok=token() oper=OP_RDSEXPR else if (tok == M_COMMA) then call s_save(OP_RDUNQUOTE, nil, nil) tok=token() oper=OP_RDSEXPR else if (tok == M_ATMARK) then call s_save(OP_RDUQTSP, nil, nil) tok=token() oper=OP_RDSEXPR else if (tok == M_ATOM) then call readstr(str) call s_return(mk_atom(trim(str))) else if (tok == M_DQUOTE) then call readstrexp(str, strlen) call s_return(mk_string(str(1:strlen))) else if (tok == M_SHARP) then call readstr(str) x=mk_const(str) if (x == nil) then call error0('Undefined sharp expression') else call s_return(x) end if else call error0('syntax error -- illegal token') end if else if (op == OP_RDLIST) then scm_args = cons(value, scm_args) tok = token() if (tok == M_COMMENT) then do call inchar(ch) if (currentline > eol) exit end do tok = token() end if if (tok == M_RPAREN) then call s_return(non_alloc_rev(nil, scm_args) ) else if (tok == M_DOT) then call s_save(OP_RDDOT, scm_args, nil) tok=token() oper=OP_RDSEXPR else call s_save(OP_RDLIST, scm_args, nil) oper=OP_RDSEXPR end if else if (op == OP_RDDOT) then if (token() /= M_RPAREN) then call error0('syntax error -- illegal dot expression') end if call s_return(non_alloc_rev(value, scm_args)) else if (op == OP_RDQUOTE) then call s_return(cons(quote, cons(value, nil))) else if (op == OP_RDQQUOTE) then call s_return(cons(qquote, cons(value, nil))) else if (op == OP_RDUNQUOTE) then call s_return(cons(unquote, cons(value, nil))) else if (op == OP_RDUQTSP) then call s_return(cons(unquotesp, cons(value, nil))) else if (op == OP_P0LIST) then if (.not.ispair(scm_args)) then if (plevel >= 0) call printatom(scm_args) call s_return(t) else if (car(scm_args) == quote .and. ok_abbrev(cdr(scm_args))) then if (plevel >= 0) write(outstr, '(a)', advance='no') '"' scm_args=cadr(scm_args) oper=OP_P0LIST else if (car(scm_args) == qquote .and. ok_abbrev(cdr(scm_args))) then if (plevel >= 0) write(outstr, '(a)', advance='no') '`' scm_args=cadr(scm_args) oper=OP_P0LIST else if (car(scm_args) == unquote .and. ok_abbrev(cdr(scm_args))) then if (plevel >= 0) write(outstr, '(a)', advance='no') ',' scm_args=cadr(scm_args) oper=OP_P0LIST else if (car(scm_args) == unquotesp .and. ok_abbrev(cdr(scm_args))) then if (plevel >= 0) write(outstr, '(a)', advance='no') ',@' scm_args=cadr(scm_args) oper=OP_P0LIST else if (plevel >= 0) write(outstr, '(a)', advance='no') '(' call s_save(OP_P1LIST, cdr(scm_args), nil) scm_args=car(scm_args) oper=OP_P0LIST end if else if (op == OP_P1LIST) then if (ispair(scm_args)) then call s_save(OP_P1LIST, cdr(scm_args), nil) if (plevel >= 0) write(outstr, '(a)', advance='no') ' ' scm_args=car(scm_args) oper=OP_P0LIST else if (scm_args /= nil) then if (plevel >= 0) then write(outstr, '(a)', advance='no') ' . ' call printatom(scm_args) end if end if if (plevel >= 0) write(outstr, '(a)', advance='no') ')' call s_return(t) end if else write(outstr, '(2a)') oper, 'is an illegal operator!' end if end subroutine opexe5 ! subroutine opexe6(op) integer :: op integer :: w, x, y ! list-length if (op == OP_LIST_LENGTH) then w=get_listlen(car(scm_args)) if (w < 0) then call error1('Not a list:', car(scm_args)) else call s_return(mk_number(w)) end if ! assq else if (op == OP_ASSQ) then x=car(scm_args) y=cadr(scm_args) do while (ispair(y)) if (.not.ispair(car(y))) then call error0('Unable to handle non-pair element') return end if if (x == caar(y)) exit y=cdr(y) end do if (ispair(y)) then call s_return(car(y)) else call s_return(f) end if ! get-closure-code else if (op == OP_GET_CLOSURE) then scm_args=car(scm_args) if (scm_args == nil) then call s_return(f) else if (isclosure(scm_args)) then call s_return(cons(lambda, car(value))) else call s_return(f) end if ! closure? else if (op == OP_CLOSUREP) then if (car(scm_args) == nil) then call s_return(f) end if call s_retbool(isclosure(car(scm_args))) else write(outstr, '(2a)') oper, 'is an illegal operator!' end if end subroutine opexe6 ! ! Mathematical functions ! subroutine opexe7(op) integer :: op integer :: x, v double precision :: rv ! sqrt if (op == OP_SQRT) then x = car(scm_args) v = int(sqrt(dfloat(get_ivalue(x)))) call s_return(mk_number(v)) ! abs else if (op == OP_ABS) then x = car(scm_args) v=abs(get_ivalue(x)) call s_return(mk_number(v)) ! exponent else if (op == OP_EXPT) then x = car(scm_args) v = cadr(scm_args) if (isnumber(x) .and. isnumber(v) .and. get_ivalue(v) >= 0) then v = get_ivalue(x) ** get_ivalue(v) call s_return(mk_number(v)) end if end if end subroutine opexe7 ! ! A few other library functions eg min, max ! subroutine opexe8(op) integer :: op integer :: x, v ! functions integer :: irandom ! min if (op == OP_MIN) then if (isnumber(car(scm_args))) then v=get_ivalue(car(scm_args)) x = cdr(scm_args) do while (x /= nil) v=min(v,get_ivalue(car(x))) x=cdr(x) end do call s_return(mk_number(v)) else call error0('Min needs at least one argument!') end if ! max else if (op == OP_MAX) then if (isnumber(car(scm_args))) then v=get_ivalue(car(scm_args)) x = cdr(scm_args) do while (x /= nil) v=max(v,get_ivalue(car(x))) x=cdr(x) end do call s_return(mk_number(v)) else call error0('Max needs at least one argument!') end if ! random else if (op == OP_RANDOM) then if (isnumber(car(scm_args))) then v=get_ivalue(car(scm_args)) call s_return(mk_number(irandom(1, v))) else call error0('Random needs an integer argument!') end if end if end subroutine opexe8 ! ! string functions ! subroutine opexe9(op) integer :: op integer :: el, i, ioerr, j, l, v, x logical :: comp, inword character (len=1) :: ch, sep character (len=20) :: str ch=' ' ! make-string if (op == OP_MKSTRING) then l=get_ivalue(car(scm_args)) if (cdr(scm_args) /= nil) then ch=get_substr(cadr(scm_args), 0, 1) end if call s_return(mk_string(repeat(ch,l))) ! string-length else if (op == OP_STRLEN) then call s_return(mk_number(get_strlen(car(scm_args)))) ! string-set! else if (op == OP_STRSET) then i=get_ivalue(cadr(scm_args)) if (i > get_strlen(car(scm_args))) then call error1('ERROR: string-set! out of bounds:', cadr(scm_args)) return end if call set_substring(car(scm_args), i, i, get_substr(caddr(scm_args), 0, 1)) call s_return(car(scm_args)) ! substring else if (op == OP_SUBSTR) then l=get_strlen(car(scm_args)) i=get_ivalue(cadr(scm_args)) if (i > l) then call error1('substring start out of bounds:', cadr(scm_args)) return end if if (cddr(scm_args) /= nil) then j=get_ivalue(caddr(scm_args)) if (j > l) then call error1('substring end out of bounds:', caddr(scm_args)) return end if else j=l end if call s_return(mk_string(get_substr(car(scm_args), i, j))) ! string-append else if (op == OP_STRAPPEND) then l=0 x=scm_args do while (x /= nil) l=l+get_strlen(car(x)) x=cdr(x) end do v=mk_string(repeat(' ',l)) i=0 l=0 x=scm_args do while (x /= nil) l=l+get_strlen(car(x)) call set_substring(v, i+1, l, get_string(car(x))) i=l x=cdr(x) end do call s_return(v) ! string-split else if (op == OP_STRSPLIT) then v=car(scm_args) l=get_strlen(v) sep=' ' inword=.false. if (cdr(scm_args) /= nil) then sep=get_string(cadr(scm_args)) inword=.true. end if x=nil j=1 do i=1, l ch=get_substr(v, i-1, i) if (ch == sep) then if (inword) then x=cons(mk_string(get_substr(v, j-1, i-1)), x) if (sep /= ' ') then j=i+1 else inword=.false. end if end if else if (.not.inword) then inword=.true. j=i end if end if end do if (inword) x=cons(mk_string(get_substr(v, j-1, l)), x) call s_return(reverse(x)) ! string=?, string?, string<=?, string>=?, substring? else if (op >= OP_STREQ .and. op <= OP_STRGE) then x = cdr(scm_args) if (.not.isstring(car(scm_args))) then call error1('Argument is not a string: ', car(scm_args)) return else str = get_string(car(scm_args)) end if do while (x /= nil) if (.not.isstring(car(x))) then call error1('Argument is not a string: ', car(x)) return else if (op == OP_STREQ) then comp=(str==get_string(car(x))) else if (op == OP_STRLT) then comp=(strget_string(car(x))) else if (op == OP_STRLE) then comp=(str<=get_string(car(x))) else if (op == OP_STRGE) then comp=(str>=get_string(car(x))) end if if (.not.comp) exit end if x=cdr(x) end do call s_retbool(comp) ! substring? else if (op == OP_STRFIND) then if (.not.isstring(car(scm_args))) then call error1('Argument is not a string: ', car(scm_args)) return else if (.not.isstring(cadr(scm_args))) then call error1('Argument is not a string: ', cadr(scm_args)) return end if i=index(get_string(cadr(scm_args)), get_string(car(scm_args))) if (i == 0) then call s_return(f) else call s_return(mk_number(i-1)) end if ! string->number else if (op == OP_STR2NUM) then str=get_string(car(scm_args)) read(str,'(i20)', iostat=ioerr) i if (ioerr/=0) i=0 call s_return(mk_number(i)) ! number->string else if (op == OP_NUM2STR) then write(str,'(i20)') get_ivalue(car(scm_args)) call s_return(mk_string(trim(adjustl(str)))) ! string->symbol else if (op == OP_STR2SYM) then call s_return(mk_symbol(get_string(car(scm_args)))) ! symbol->string else if (op == OP_SYM2STR) then call s_return(mk_string(get_string(caar(scm_args)))) end if end subroutine opexe9 ! ! Nonstandard library additions for system interface such as system, readline ! subroutine opexe10(op) use locus_data use locus_types integer :: op integer, parameter :: MISS = -9999 integer :: i, j, l, tmp, typ, v, x character (len=20) :: str character (len=256) :: buff ! functions logical :: strfind ! system command from within scheme str=' ' if (op == OP_SYSTEM) then if (.not.isstring(car(scm_args))) then call error0('Error -- first argument must be a string') call s_return(f) else call system(get_string(car(scm_args))) call s_return(t) end if else if (op == OP_MKPORT) then if (isstring(car(scm_args))) then buff=get_string(car(scm_args)) end if write(outstr, '(3a)') 'Opening "', trim(buff), '"' if (nports < MAXPORT) then nports=nports+1 iport=portaddress(nports) open(iport, file=trim(buff), status='unknown') call s_return(mk_port(iport, trim(buff))) else call error0('Too many open files!') call s_return(f) end if else if (op == OP_CLPORT) then if (isport(car(scm_args))) then iport=get_ivalue(car(scm_args)) j=test_port(iport) if (j > 0) then close(iport) tmp=portaddress(nports) portaddress(nports)=portaddress(j) portaddress(j)=tmp nports=nports-1 call set_ivalue(car(scm_args), 0) call s_return(t) else call error0('Closed port!') call s_return(f) end if else call error0('Not a port!') call s_return(f) end if else if (op == OP_RDLINE) then if (car(scm_args) == nil) then read(*, '(a)', iostat=ioerr) buff if (ioerr==0) then call s_return(mk_string(trim(buff))) else call s_return(f) end if else if (isport(car(scm_args))) then iport=get_ivalue(car(scm_args)) j=test_port(iport) if (j > 0) then read(iport, '(a)', iostat=ioerr) buff if (ioerr==0) then call s_return(mk_string(trim(buff))) else call s_return(f) end if else call error0('Closed port!') call s_return(f) end if else call error0('Not a port!') call s_return(f) end if ! run a Sib-pair command else if (op == OP_RUNCMD) then l=0 x=scm_args do while (x /= nil) if (isnumber(car(x))) then write(str, '(i20)') get_ivalue(car(x)) str=adjustl(str) l=l+len_trim(str)+1 else l=l+get_strlen(car(x))+1 end if x=cdr(x) end do v=mk_string(repeat(' ',l)) i=0 l=0 x=scm_args do while (x /= nil) if (isnumber(car(x))) then write(str, '(i20)') get_ivalue(car(x)) str=adjustl(str) l=l+len_trim(str) call set_substring(v, i+1, l, trim(str)) else l=l+get_strlen(car(x)) call set_substring(v, i+1, l, get_string(car(x))) end if l=l+1 i=l call set_substring(v, i, i, ';') x=cdr(x) end do commands=get_string(v) // trim(commands) if (ilevel /= 0) ilevold=ilevel ilevel=0 call s_return(t) ! Scheme specific locus list else if (op == OP_LSLOCI) then typ=0 if (isnumber(car(scm_args))) then typ = get_ivalue(car(scm_args)) else if (isstring(car(scm_args)) .or. issymbol(car(scm_args))) then if (isstring(car(scm_args))) then x=car(scm_args) else x=caar(scm_args) end if typ=loccode(get_substr(x, 0, 1)) if (typ==LOC_DEL .and. len_trim(get_string(x)) == 2) then i=loccode(get_substr(x, 1, 2)) if (i > 0 .and. i < LOC_DEL) typ=typ+i end if end if x=nil if (typ == 0) then do i=1, nloci x=cons(mk_string(trim(loc(i))), x) end do else if (typ == 10) then do i=1, nloci if (loctyp(i) > typ) then x=cons(mk_string(trim(loc(i))), x) end if end do else do i=1, nloci if (loctyp(i) == typ) then x=cons(mk_string(trim(loc(i))), x) end if end do end if call s_return(reverse(x)) ! number of loci else if (op == OP_NLOCI) then call s_return(mk_number(nloci)) ! Return information about locus else if (op >= OP_LOCNAM .and. op <= OP_LOCRANK) then i=0 if (isnumber(car(scm_args))) then i = get_ivalue(car(scm_args)) else if (isstring(car(scm_args))) then str=trim(get_string(car(scm_args))) do j=1, nloci if (str == loc(j)) then i=j exit end if end do end if if (i > 0 .and. i <= nloci) then if (op == OP_LOCNAM) then call s_return(mk_string(trim(loc(i)))) else if (op == OP_LOCTYP) then i=loctyp(i) if (i < LOC_DEL) then str=typloc(i) else str='d' str(2:2)=typloc(i-LOC_DEL) end if call s_return(mk_string(trim(str))) else if (op == OP_LOCORD) then call s_return(mk_number(i)) else if (op == OP_LOCNOT) then call s_return(mk_string(trim(locnotes(i)))) else if (op == OP_LOCRANK) then v=0 j=1 do while (wloc(j) > 0) if (wloc(j)==i) then v=j exit end if j=j+1 end do call s_return(mk_number(v)) end if else call s_return(f) end if ! Scheme specific help else if (op == OP_APROPOS) then str=' ' if (isstring(car(scm_args))) then str='*' // trim(get_string(car(scm_args))) // '*' end if tmp=oblist do while (tmp /= nil) if (strfind(trim(str), get_string(caar(tmp)), 1)) then write(outstr, '(a)') get_string(caar(tmp)) end if tmp=cdr(tmp) end do call s_return(t) else if (op == OP_HELP) then write(outstr, '(a)') & 'Sib-pair Scheme is a minimal scheme interpreter.', & 'It implements integer arithmetic and strings.', & '(apropos ) lists commands containing that string.', & '(read-line []) reads in next line from stdin or open file.', & '(string-split []) splits string on white space or optional char.', & '(substring? ) returns start of substring in string.', & '(ls []) creates a list of locus names (of given type "adhmqx").', & '(nloci) returns total number of loci.', & '(loc ) returns locus at that position in the locus list.', & '(locord ) returns position of a locus in the locus list.', & '(locnotes ) returns notes for a locus.', & '(loctyp ) evaluates type of a locus ("adhmqx").', & '(run ...) stores Sib-pair commands to the buffer', & 'for evaluation once you return to the usual Sib-pair prompt.' call s_return(t) end if end subroutine opexe10 ! ! Initialization of internal keywords ! subroutine mk_syntax(op, nam) integer :: op character (len=*) :: nam integer :: x x=cons(mk_string(nam), nil) call set_type(x, ior(T_SYNTAX, T_SYMBOL)) call set_syntaxnum(x, op) oblist=cons(x, oblist) end subroutine mk_syntax ! subroutine mk_proc(op, nam) integer :: op character (len=*) :: nam integer :: x, y x=mk_symbol(nam) y=getcell(nil, nil) call set_type(y, ior(T_PROC, T_ATOM)) call set_ivalue(y, op) call set_car(global_env, cons(cons(x, y), car(global_env))) end subroutine mk_proc ! ! Initiate global environment ! subroutine init_vars_global() integer :: x global_env=cons(nil, nil) x=mk_symbol('else') call set_car(global_env, cons(cons(x, t), car(global_env))) end subroutine init_vars_global ! ! Initiate syntax ! subroutine init_syntax() call mk_syntax(OP_LAMBDA, 'lambda') call mk_syntax(OP_QUOTE, 'quote') call mk_syntax(OP_DEF0, 'define') call mk_syntax(OP_IF0, 'if') call mk_syntax(OP_BEGIN, 'begin') call mk_syntax(OP_SET0, 'set!') call mk_syntax(OP_LET0, 'let') call mk_syntax(OP_LET0AST, 'let*') call mk_syntax(OP_LET0REC, 'letrec') call mk_syntax(OP_COND0, 'cond') call mk_syntax(OP_DELAY, 'delay') call mk_syntax(OP_AND0, 'and') call mk_syntax(OP_OR0, 'or') call mk_syntax(OP_C0STREAM, 'cons-stream') call mk_syntax(OP_0MACRO, 'macro') call mk_syntax(OP_CASE0, 'case') end subroutine init_syntax ! ! Initiate procedures ! subroutine init_procs() call mk_proc(OP_PEVAL, 'eval') call mk_proc(OP_PAPPLY, 'apply') call mk_proc(OP_CONTINUATION, 'call-with-current-continuation') call mk_proc(OP_FORCE, 'force') call mk_proc(OP_CAR, 'car') call mk_proc(OP_CDR, 'cdr') call mk_proc(OP_CONS, 'cons') call mk_proc(OP_SETCAR, 'set-car!') call mk_proc(OP_SETCDR, 'set-cdr!') call mk_proc(OP_ADD, '+') call mk_proc(OP_SUB, '-') call mk_proc(OP_MUL, '*') call mk_proc(OP_DIV, '/') call mk_proc(OP_INTDIV, 'quotient') call mk_proc(OP_REM, 'remainder') call mk_proc(OP_MOD, 'modulo') call mk_proc(OP_NOT, 'not') call mk_proc(OP_BOOL, 'boolean?') call mk_proc(OP_ISINT, 'integer?') call mk_proc(OP_SYMBOL, 'symbol?') call mk_proc(OP_NUMBER, 'number?') call mk_proc(OP_STRING, 'string?') call mk_proc(OP_PROC, 'procedure?') call mk_proc(OP_PAIR, 'pair?') call mk_proc(OP_EQV, 'eqv?') call mk_proc(OP_EQ, 'eq?') call mk_proc(OP_NULL, 'null?') call mk_proc(OP_ZEROP, 'zero?') call mk_proc(OP_POSP, 'positive?') call mk_proc(OP_NEGP, 'negative?') call mk_proc(OP_NUMEQ, '=') call mk_proc(OP_LESS, '<') call mk_proc(OP_GRE, '>') call mk_proc(OP_LEQ, '<=') call mk_proc(OP_GEQ, '>=') call mk_proc(OP_READ, 'read') call mk_proc(OP_WRITE, 'write') call mk_proc(OP_DISPLAY, 'display') call mk_proc(OP_NEWLINE, 'newline') call mk_proc(OP_ERR0, 'error') call mk_proc(OP_REVERSE, 'reverse') call mk_proc(OP_APPEND, 'append') call mk_proc(OP_GC, 'gc') call mk_proc(OP_GCVERB, 'gc-verbose') call mk_proc(OP_NEWSEGMENT, 'new-segment') call mk_proc(OP_LIST_LENGTH, 'length') call mk_proc(OP_ASSQ, 'assq') call mk_proc(OP_GET_CLOSURE, 'get-closure-code') call mk_proc(OP_CLOSUREP, 'closure?') call mk_proc(OP_QUIT, 'quit') call mk_proc(OP_SQRT, 'sqrt') call mk_proc(OP_ABS, 'abs') call mk_proc(OP_EXPT, 'expt') call mk_proc(OP_MIN, 'min') call mk_proc(OP_MAX, 'max') call mk_proc(OP_RANDOM, 'random') call mk_proc(OP_MKSTRING, 'make-string') call mk_proc(OP_STRLEN, 'string-length') call mk_proc(OP_STRSET, 'string-set!') call mk_proc(OP_SUBSTR, 'substring') call mk_proc(OP_STRAPPEND, 'string-append') call mk_proc(OP_STRSPLIT, 'string-split') call mk_proc(OP_STREQ, 'string=?') call mk_proc(OP_STRLT, 'string?') call mk_proc(OP_STRLE, 'string<=?') call mk_proc(OP_STRGE, 'string>=?') call mk_proc(OP_STRFIND, 'substring?') call mk_proc(OP_STR2NUM, 'string->number') call mk_proc(OP_NUM2STR, 'number->string') call mk_proc(OP_STR2SYM, 'string->symbol') call mk_proc(OP_SYM2STR, 'symbol->string') call mk_proc(OP_SYSTEM, 'system') call mk_proc(OP_MKPORT, 'open-input-file') call mk_proc(OP_CLPORT, 'close-input-port') call mk_proc(OP_RDLINE, 'read-line') call mk_proc(OP_RUNCMD, 'run') call mk_proc(OP_LSLOCI, 'ls') call mk_proc(OP_NLOCI, 'nloci') call mk_proc(OP_LOCNAM, 'loc') call mk_proc(OP_LOCTYP, 'loctyp') call mk_proc(OP_LOCORD, 'locord') call mk_proc(OP_LOCNOT, 'locnotes') call mk_proc(OP_LOCRANK, 'locrank') call mk_proc(OP_APROPOS, 'apropos') call mk_proc(OP_HELP, 'help') end subroutine init_procs ! ! Inlined init.scm ! subroutine init_scm() lin='(define nil #f) (define t #t)' call repl_scheme(3,0) lin='(define (caar x) (car (car x)))' call repl_scheme(3,0) lin='(define (cadr x) (car (cdr x)))' call repl_scheme(3,0) lin='(define (cdar x) (cdr (car x)))' call repl_scheme(3,0) lin='(define (cddr x) (cdr (cdr x)))' call repl_scheme(3,0) lin='(define (caaar x) (car (car (car x))))' call repl_scheme(3,0) lin='(define (caadr x) (car (car (cdr x))))' call repl_scheme(3,0) lin='(define (cadar x) (car (cdr (car x))))' call repl_scheme(3,0) lin='(define (caddr x) (car (cdr (cdr x))))' call repl_scheme(3,0) lin='(define (cdaar x) (cdr (car (car x))))' call repl_scheme(3,0) lin='(define (cdadr x) (cdr (car (cdr x))))' call repl_scheme(3,0) lin='(define (cddar x) (cdr (cdr (car x))))' call repl_scheme(3,0) lin='(define (cdddr x) (cdr (cdr (cdr x))))' call repl_scheme(3,0) lin='(define call/cc call-with-current-continuation)' call repl_scheme(3,0) lin='(define (list . x) x)' call repl_scheme(3,0) ! lin='(define (map proc list) (if (pair? list) (cons (proc (car list)) (map proc (cdr list)))))' lin='(define (unzip1-with-cdr . lists) ' // & '(unzip1-with-cdr-iterative lists ' // char(39) // '()' // & char(39) // '()))' call repl_scheme(3,0) lin='(define (unzip1-with-cdr-iterative lists cars cdrs) ' // & '(if (null? lists) (cons cars cdrs) ' // & '(let ((car1 (caar lists)) (cdr1 (cdar lists)))' // & '(unzip1-with-cdr-iterative (cdr lists) ' // & ' (append cars (list car1)) (append cdrs (list cdr1))))))' call repl_scheme(3,0) lin='(define (map proc . lists) (if (null? lists) (apply proc)' // & '(if (null? (car lists)) ' // char(39) // '() ' // & '(let* ((unz (apply unzip1-with-cdr lists)) (cars (car unz)) ' // & '(cdrs (cdr unz))) (cons (apply proc cars) ' // & '(apply map (cons proc cdrs)))))))' call repl_scheme(3,0) lin='(define (for-each proc list) (if (pair? list) (begin (proc (car list)) (for-each proc (cdr list))) #t ))' call repl_scheme(3,0) lin='(define (list-tail x k) (if (zero? k) x (list-tail (cdr x) (- k 1))))' call repl_scheme(3,0) lin='(define (list-ref x k) (if (null? x) x (car (list-tail x k))))' call repl_scheme(3,0) lin='(define (atom? x) (not (pair? x)))' call repl_scheme(3,0) lin='(define (memq obj lst) (cond ((null? lst) #f) ((eq? obj (car lst)) lst) (else (memq obj (cdr lst)))))' call repl_scheme(3,0) lin='(define (equal? x y) (if (pair? x) (and (pair? y)' // & '(equal? (car x) (car y)) (equal? (cdr x) (cdr y)))' // & '(and (not (pair? y)) (eqv? x y))))' call repl_scheme(3,0) lin='(define (even? x) (if (integer? x) (zero? (remainder x 2))))' call repl_scheme(3,0) lin='(define (odd? x) (if (integer? x) (= (remainder x 2) 1)))' call repl_scheme(3,0) lin='(define (gcd a b) (let ((aa (abs a)) (bb (abs b)))' // & '(if (= bb 0) aa (gcd bb (remainder aa bb)))))' call repl_scheme(3,0) lin='(define (lcm a b) (if (or (= a 0) (= b 0)) 0 ' // & '(abs (* (quotient a (gcd a b)) b))))' call repl_scheme(3,0) lin='(define (generic-member cmp obj lst) (cond ((null? lst) #f)' // & '((cmp obj (car lst)) lst)' // & '(else (generic-member cmp obj (cdr lst)))))' call repl_scheme(3,0) lin='(define (memq obj lst) (generic-member eq? obj lst))' call repl_scheme(3,0) lin='(define (memv obj lst) (generic-member eqv? obj lst))' call repl_scheme(3,0) lin='(define (member obj lst) (generic-member equal? obj lst))' call repl_scheme(3,0) lin='(define (generic-assoc cmp obj alst)' // & '(cond ((null? alst) #f) ((cmp obj (caar alst))' // & '(car alst)) (else (generic-assoc cmp obj (cdr alst)))))' call repl_scheme(3,0) lin='(define (assq obj alst) (generic-assoc eq? obj alst))' call repl_scheme(3,0) lin='(define (assv obj alst) (generic-assoc eqv? obj alst))' call repl_scheme(3,0) lin='(define (assoc obj alst) (generic-assoc equal? obj alst))' call repl_scheme(3,0) lin='(define (exit) ((run "exit") (quit)))' call repl_scheme(3,0) end subroutine init_scm ! ! Initiate procedures ! subroutine init_globals() call init_vars_global() call init_syntax() call init_procs() ! init global pointers to special symbols lambda=mk_symbol('lambda') quote=mk_symbol('quote') qquote=mk_symbol('quasiquote') unquote=mk_symbol('unquote') unquotesp=mk_symbol('unquote-splicing') call init_scm() end subroutine init_globals ! ! Error handling ! subroutine error0(s) character (len=*) :: s scm_args=cons(mk_string(s), nil) oper=OP_ERR0 write(*,'(/5a)') ' At: "', lin(1:currentline), '^', & trim(lin(currentline+1:eol)) , '"' write(*,'(8x,a)') s end subroutine error0 ! subroutine error1(s, a) integer :: a character (len=*) :: s scm_args=cons(a, nil) scm_args=cons(mk_string(s), scm_args) oper=OP_ERR0 write(*,'(/5a)') ' At: "', lin(1:currentline), '^', & trim(lin(currentline+1:eol)) , '"' end subroutine error1 ! subroutine init_scheme() call setup_mem(100) call init_globals() end subroutine init_scheme ! ! Scheme read-eval-print loop ! subroutine repl_scheme(inline, ple) integer, intent(in) :: inline integer, intent(in) :: ple integer :: i, it, mlevel, op, plevel plevel=0 mlevel=inline oper = OP_T0LVL prompt_string='%% ' currentline=0 if (mlevel==1) then lin=' ' else if (mlevel == 2) then i=1 do while (lin(i:i) == ' ') i=i+1 end do do while (lin(i:i) /= ' ') i=i+1 end do i=i+1 lin=lin(i:len_trim(lin)) // ' (quit)' if (ple < -1) plevel=-1 else plevel=-1 lin=trim(lin) // ' (quit)' end if eol=len_trim(lin) ! do it=it+1 op=oper if (op == OP_ERR0 .or. op == OP_ERR1) then if (mlevel==2) exit mlevel=1 end if if (op == OP_QUIT .or. irupt > 0) then if (inline==1) write(*, '(a)') 'Leaving scheme!' exit else if (op >= OP_LOAD .and. op <= OP_LET2AST) then call opexe0(op, inline) else if (op >= OP_LET0REC .and. op <= OP_CONTINUATION) then call opexe1(op) else if (op >= OP_ADD .and. op <= OP_SETCDR) then call opexe2(op) else if (op >= OP_NOT .and. op <= OP_EQV) then call opexe3(op) else if (op >= OP_FORCE .and. op <= OP_NEWSEGMENT) then call opexe4(op, plevel) else if (op >= OP_RDSEXPR .and. op <= OP_P1LIST) then call opexe5(op, plevel) else if (op >= OP_LIST_LENGTH .and. op <= OP_MACROP) then call opexe6(op) else if (op >= OP_EXP .and. op <= OP_EXPT) then call opexe7(op) else if (op >= OP_MIN .and. op <= OP_RANDOM) then call opexe8(op) else if (op >= OP_MKSTRING .and. op <= OP_STR2SYM) then call opexe9(op) else if (op >= OP_SYSTEM .and. op <= OP_HELP) then call opexe10(op) else write(*, '(a)') 'Bad op code! Exiting!' exit end if end do prompt_string='>> ' end subroutine repl_scheme end module scheme_lang ! ! Simple regression formula structure and parser ! ! formula and design matrix ! formula is: a b c a*b a*c b*c a*b*c ! T1 1 2 3 1 1 2 1 ! T2 . . . 2 3 3 2 ! T3 . . . . . . 3 ! TERMDIM 1 1 1 2 2 2 3 ! ! Effects 1 2 3 ! NLEV n1 n2 n3 ! STA 1 n1+1 n1+n2+1 ! FIN n1 n1+n2 n1+n2+n3 ! INFORM 1 1 1 ! module formula_class type formula_data logical :: intercept ! add an intercept term integer :: nterms ! number of model terms integer :: maxlev ! maximum number of effects in an interaction integer :: neff ! number of main effects integer :: mainlen ! number of main effects columns integer :: maxrows ! maximum rows in design matrix integer, dimension(:), allocatable :: effects ! included effects integer, dimension(:), allocatable :: nlev ! no. levels of each effect integer, dimension(:), allocatable :: inform ! if/how included in model formula integer, dimension(:), allocatable :: sta ! start of effect in main effects design matrix integer, dimension(:), allocatable :: fin ! end of effect in main effects design matrix integer, dimension(:,:), allocatable :: termlist ! matrix storing term components integer, dimension(:), allocatable :: termdim ! interaction order integer, dimension(:), allocatable :: termlev ! no. levels of each interaction term integer :: designcols ! total number of design matrix columns end type formula_data contains ! ! read commands and write appropriate formula structure ! subroutine create_form(sta, fin, terms, nloci, loc, formula) integer :: sta, fin character (len=*), dimension(:), intent(in) :: terms integer, intent(in) :: nloci character (len=20), dimension(:), intent(in) :: loc type (formula_data) :: formula integer, dimension(fin-sta+1) :: tmpeff integer :: fpos, i, lev, locnum, pos, prev logical :: fnd, interact character (len=1) :: cross = '*', minus = '-', plus = '+' ! functions integer isinenv formula%intercept=.true. formula%nterms=0 formula%maxrows=0 formula%maxlev=1 formula%neff=0 interact=.false. lev=1 prev=0 pos=sta ! prev=0 between terms ! prev=1 reading an interaction (prev symbol = *) ! prev=2 read a variable name do while (pos <= fin) if (terms(pos)==cross) then if (prev == 2) then lev=lev+1 interact=.true. if (lev > formula%maxlev) formula%maxlev=lev prev=1 else write(*,'(3a)') 'ERROR: In model formula "*" requires two arguments.' end if else if (terms(pos)==plus) then prev=0 else if (terms(pos)==minus) then if (pos < fin) then if (terms(pos+1) == '1') then pos=pos+1 else write(*,'(3a)') 'ERROR: In model formula "-', trim(terms(pos+1)), '" not allowed.' end if end if prev=0 else locnum=isinenv(terms(pos)(1:20), nloci, loc) if (locnum /= 0) then fnd=.false. do i=1, formula%neff if (locnum == tmpeff(i)) then fnd=.true. exit end if end do if (.not.fnd) then formula%neff=formula%neff+1 tmpeff(formula%neff)=locnum end if if (interact) then prev=2 interact=.false. else formula%nterms=formula%nterms+1 prev=2 lev=1 end if else write(*,'(3a)') 'ERROR: In model formula "', trim(terms(pos)), '" is not a locus.' end if end if pos=pos+1 end do allocate(formula%effects(formula%neff)) allocate(formula%nlev(formula%neff)) allocate(formula%inform(formula%neff)) allocate(formula%sta(formula%neff)) allocate(formula%fin(formula%neff)) allocate(formula%termlist(formula%nterms, formula%maxlev)) allocate(formula%termdim(formula%nterms)) allocate(formula%termlev(formula%nterms)) formula%effects=tmpeff(1:formula%neff) formula%mainlen=0 formula%nlev=0 formula%inform=0 formula%sta=0 formula%fin=0 formula%termlist=0 formula%termdim=0 formula%termlev=0 fpos=0 interact=.false. prev=0 pos=sta do while (pos <= fin) if (terms(pos)==cross) then if (prev == 2) then formula%termdim(fpos)=formula%termdim(fpos)+1 interact=.true. prev=1 end if else if (terms(pos)==plus) then prev=0 else if (terms(pos)==minus) then if (pos < fin) then if (terms(pos+1) == '1') then formula%intercept=.false. pos=pos+1 end if end if prev=0 else locnum=isinenv(terms(pos)(1:20), nloci, loc) if (locnum /= 0) then if (prev /= 1) then prev=2 fpos=fpos+1 formula%termdim(fpos)=1 end if do i=1, formula%neff if (locnum == formula%effects(i)) then formula%termlist(fpos,formula%termdim(fpos))=i exit end if end do prev=2 if (interact) then interact=.false. end if end if end if pos=pos+1 end do do fpos=1, formula%nterms if (formula%termdim(fpos)==1) then formula%inform(formula%termlist(fpos,1))=1 end if end do end subroutine create_form ! subroutine cleanup_form(formula) type (formula_data) :: formula formula%intercept=.false. formula%nterms=0 formula%maxlev=0 formula%neff=0 formula%mainlen=0 formula%maxrows=0 formula%designcols=0 deallocate(formula%effects, formula%nlev, formula%sta, formula%fin) deallocate(formula%termlist, formula%termdim, formula%termlev) end subroutine cleanup_form ! ! If levels available, evaluate total number of model parameters ! subroutine sumcols_form(formula) type (formula_data) :: formula integer :: i, j, df, hasmain, lpos, off, pos logical :: anymain off=0 if (formula%intercept) off=-1 formula%designcols=formula%designcols-off formula%maxrows=1 pos=0 do i=1, formula%neff formula%sta(i)=pos+1 formula%fin(i)=pos+formula%nlev(i) formula%mainlen=formula%mainlen+formula%nlev(i) formula%maxrows=formula%maxrows*formula%nlev(i) pos=pos+formula%nlev(i) end do do i=1, formula%nterms anymain=.false. df=1 do j=1, formula%termdim(i) lpos=formula%termlist(i,j) if (formula%nlev(lpos) > 0) then hasmain=0 if (formula%inform(lpos)==1) then anymain=.true. hasmain=1 end if df=df*(formula%nlev(lpos)-hasmain) end if end do if (.not.anymain) df=df+off formula%termlev(i)=df formula%designcols=formula%designcols+df end do end subroutine sumcols_form ! subroutine show_form(formula) type (formula_data) :: formula integer :: i write(*,*) 'Intercept : ', formula%intercept write(*,*) 'Nterms : ', formula%nterms write(*,*) 'Maxlev : ', formula%maxlev write(*,*) 'Neff : ', formula%neff write(*,*) 'Mainlen : ', formula%mainlen write(*,*) 'Maxrows : ', formula%maxrows write(*,*) 'DesignCols: ', formula%designcols write(*,*) 'Effects : ', formula%effects write(*,*) 'Nlev : ', formula%nlev write(*,*) 'Sta : ', formula%sta write(*,*) 'Fin : ', formula%fin write(*,*) 'InFormula : ', formula%inform write(*,*) 'Termlist : ', formula%termlist(1:formula%nterms,1) do i=2, formula%maxlev write(*,*) ' ', formula%termlist(1:formula%nterms,i) end do write(*,*) 'Termdim : ', formula%termdim write(*,*) 'Termlev : ', formula%termlev end subroutine show_form end module formula_class ! ! Mixed model analysis parameters ! ! linkf=link function 1=identity 2=logit 3=probit 4=MFT 5=log ! modtyp=likelihood family (1=gaussian, 2=binom