Seg fault with "wrapped" code

From: Ángel G. Muñoz <agmunoz_at_nyahnyahspammersnyahnyah>
Date: Thu Jan 26 2012 - 21:22:34 MST

Hello,

I have successfully wrapped a Fortran code using WRAPIT, but it looks like I don't understand the way some multidimensional arrays are passed to the fortran code.

The main issue is that NCL claims that the sizes of the arrays are commuted for some arrays but not for others. It just makes no sense. But when I solve all the problems NCL reports, then I obtain a seg fault.

The script is a simple modification of rip_cape.f, for computing the normalized CAPE and CIN (Blanchard, 1998).

I am attaching the main lines of the ncl script and the fortran code.

Thanks in advance

 load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
  load "$NCARG_ROOT/lib/ncarg/nclscripts/wrf/WRFUserARW.ncl"
  external NCAPENCIN "./NCAPE_NCIN.so"

begin
  ano=2006
  archs=systemfunc ("ls Cata2006_CAPEvars.nc Cata2007_CAPEvars.nc Cata2008_CAPEvars.nc") ;Lista de archivos a procesar

  print("Se procesaran los siguientes archivos:")
  print(archs)

  tam=dimsizes(archs)

  do i=0,tam-1 ;Hacemos esto para cada archivo de la lista a procesar
    ;Archivo de salida
    archsalid = ""+ano+".nc"
    system ("rm -Rf "+archsalid)

    ;Lecturas varias
    a = addfile(archs(i),"r")
    times = wrf_user_list_times(a)
    ntimes = dimsizes(times)
    anol = stringtoint(str_get_cols(times, 0, 3))
    mesl = stringtoint(str_get_cols(times, 5, 6))

    lat2d = a->XLAT(0,:,:)
    lon2d = a->XLONG(0,:,:)
    T = a->T
    P = a->P
    PB = a->PB
    QV = a->QVAPOR
    PH = a->PH
    PHB = a->PHB
    HGT = a->HGT(0,:,:)
    PSFC = a->PSFC
    T = T + 300.
    P = P + PB
    tk = wrf_tk( P , T )
    PH = PH + PHB
    z = wrf_user_unstagger(PH,PH@stagger)
    z = z/9.81

;Transformacion de coords curvilineas a cartesianas (req. por subrutina)
    ;lat = fspan(8.,11.,109)
    ;lon = fspan(-73.,-70.,182)
    ;Pr = rcm2rgrid(lat2d,lon2d,P(0,:,:),lat,lon,0)
    ;TKr = rcm2rgrid(lat2d,lon2d,tk(0,:,:),lat,lon,0)
    ;QVr = rcm2rgrid(lat2d,lon2d,QV(0,:,:),lat,lon,0)
    ;zr = rcm2rgrid(lat2d,lon2d,z(0,:,:),lat,lon,0)
    ;HGTr = rcm2rgrid(lat2d,lon2d,HGT(0,:,:),lat,lon,0)
    ;PSFCr = rcm2rgrid(lat2d,lon2d,PSFC(0,:,:),lat,lon,0)

;Paja adicional (pero necesaria por la subrutina):
    I3DFLAG =1
    TER_FOLLOW=1
    PSAFILE ="psadilookup.dat"
    MIY =dimsizes(PSFC(0,:,0))
    MJX =dimsizes(PSFC(0,0,:))
    MKZH =27
    CAPE =new((/MIY,MJX,MKZH/),"float")
    CIN =new((/MIY,MJX,MKZH/),"float")
    NCAPE =new((/MIY,MJX,MKZH/),"float")
    NCIN =new((/MIY,MJX,MKZH/),"float")

    P!0="time"
    P!1="lev"
    P!2="lat"
    P!3="lon"
    tk!0="time"
    tk!1="lev"
    tk!2="lat"
    tk!3="lon"
    QV!0="time"
    QV!1="lev"
    QV!2="lat"
    QV!3="lon"
    z!0="time"
    z!1="lev"
    z!2="lat"
    z!3="lon"
; HGT!0="time"
    HGT!0="lat"
    HGT!1="lon"
    PSFC!0="time"
    PSFC!1="lat"
    PSFC!2="lon"

    Pf =P(time|:,lat|:,lon|:,lev|:)
    tkf =tk(time|:,lat|:,lon|:,lev|:)
    QVf =QV(time|:,lat|:,lon|:,lev|:)
    zf =z(time|:,lat|:,lon|:,lev|:)
    HGTf =HGT(lon|:,lat|:)
    PSFCf =PSFC(time|:,lon|:,lat|:)

printVarSummary(HGTf)
printVarSummary(Pf)
print(MIY)
print(MJX)

    ;Calculo NCAPE y NCIN
    print("Inicio de calculo de NCAPE y NCIN")

    do ti=0,ntimes-1
print(ti)
      Pf2=Pf(ti,:,:,:)
      tkf2=tkf(ti,:,:,:)
      QVf2=QVf(ti,:,:,:)
      zf2=zf(ti,:,:,:)
      PSFCf2=PSFCf(ti,:,:)

; NCAPENCIN::DCAPECALC3D(Pr,TKr,QVr,zr,HGTr,PSFCr,CAPE,CIN,NCAPE,NCIN,MIY,MJX,MKZH,I3DFLAG,TER_FOLLOW,PSAFILE)
      NCAPENCIN::DCAPECALC3D(Pf2,tkf2,QVf2,zf2,HGTf,PSFCf2,CAPE,CIN,NCAPE,NCIN,MIY,MJX,MKZH,I3DFLAG,TER_FOLLOW,PSAFILE)
      NCAPE4D(ti,:,:,:)=NCAPE
      NCIN4D(ti,:,:,:) =NCIN
      delete(CAPE)
      delete(NCAPE)
      delete(CIN)
      delete(NCIN)
    end do

    print("NCAPE y NCIN calculados")
    ndim = dimsizes(NCAPE4D)
    print(NCAPE)
 

Ángel G. Muñoz S.
Coordinador Eje Geociencias
Centro de Modelado Científico (CMC)
La Universidad del Zulia
VENEZUELA
http://cmc.org.ve

_______________________________________________
ncl-talk mailing list
List instructions, subscriber options, unsubscribe:
http://mailman.ucar.edu/mailman/listinfo/ncl-talk
Received on Fri Jan 27 09:22:47 2012

This archive was generated by hypermail 2.1.8 : Thu Feb 02 2012 - 03:10:31 MST