function sunrot, Data, date, time, date1, time1, LIMB, $
LATITUDE=LATITUDE, SHOW=SHOW, outside=outside, missing=missing

;+	Compensation of differential rotation of the Sun.
;	If keyword parameter LATITUDE is present, the 'solid' rotation
;	is performed, and LATITUDE is interpreted is the reference one.
;	Example:
;
;	Data=rfits(pickfile(/read,filt='*.fts'))
;	Sz=size(Data)
;	window,/free,xsize=Sz(1),ysize=Sz(2),xpos=0,ypos=0
;	LIMB=[Center_X, Center_Y, Radius]
;	Date='01 01 93'
;	Time='02:07:35'
;	Date1='02 01 93'
;	Time='05:17:47'
;	Rotated_Sun=sunrot(Data,date,time,date1,time1,LIMB)
;	tvscl,Rotated_Sun
;-


if n_params() lt 6 then message, 'Insufficient number of arguments'
if n_elements(missing) le 0 then missing=0

dt=time_difference(date, time, date1, time1)/3600d0/24

WIDGET_CONTROL,/hour

Sz=size(Data)

type=Sz(Sz(0)+1)

suneph,date,time,SUN

Lonmin=-89.
Lonmax=89.
Latmin=float(-89+SUN.B0*!Radeg)
Latmax=float(89+SUN.B0*!Radeg)

Centre=float(Limb([0,1]))
Radius=float(Limb(2))

;*****************************************************************

window,/free,xs=Sz(1),ys=Sz(2), pixmap=1-keyword_set(show)
N_win=!d.window

if keyword_set(show) then tvscl,data

	if strmid(!version.release,0,1) lt 5 then begin

!x.style=(!y.style=1)

!x.range=[-Centre(0),!d.x_size-Centre(0)]/Radius
!y.range=[-Centre(1),!d.y_size-Centre(1)]/Radius

map_set,float(SUN.B0*!Radeg),0,0, /grid, $
        /ortho,/noerase,pos=[0,0,1,1],/nobor,latdel=10,londel=10,col=100

!x.style=(!y.style=(!x.range=(!y.range=0)))

	endif else begin

map_set,float(SUN.B0*!Radeg),0,0, /ortho,/noerase,pos=[0,0,1,1],/nobor

!x.s=[Centre(0), Radius] / float(!d.x_size)
!y.s=[Centre(1), Radius] / float(!d.y_size)

P_clip_save=!P.clip
!p.clip = [0, 0, !d.x_size, !d.y_size]

MAP_GRID, latdel=10, londel=10, col=100

!P.clip=P_clip_save


	endelse

Nm=512*2
arg=(findgen(Nm)-Nm/4)*(!Pi*2./Nm)

ax=cos(arg)
ay=sin(arg)

xcir=ax*(Radius-2)+Centre(0)
ycir=ay*(Radius-2)+Centre(1)

Disk_index=(polyfillv(xcir,ycir,Sz(1),Sz(2)))

xcir=ax*Radius+Centre(0)
ycir=ay*Radius+Centre(1)

median_value=100b

z=make_array(Sz(1),Sz(2), type=type)

z(polyfillv(xcir,ycir,Sz(1),Sz(2)))=median_value

Planecoord=(subs1to2(Disk_index,dimension=[Sz(1),Sz(2)]))

if keyword_set(outside) then Disk=z eq 0

Sphcoord=(convert_coord(Planecoord,/dev,/to_data))([0,1],*)

bad_ind=where(			$
		(Sphcoord(0,*) le Lonmin) and 	$
		(Sphcoord(0,*) ge Lonmax) and 	$
		(Sphcoord(1,*) le Latmin) and 	$
		(Sphcoord(1,*) ge Latmax)	)


if bad_ind(0) ge 0 then Sphcoord(*,bad_ind)=0


if n_elements(LATITUDE) le 0 then $
	Sphcoord(0,*)=Sphcoord(0,*)-difrot(dt,Sphcoord(1,*),/days,/degree) else $
	Sphcoord(0,*)=Sphcoord(0,*)-difrot(dt,Sphcoord(1,*),/days,/degree, $
		lat=LATITUDE)

NewPlanecoord=(convert_coord(temporary(Sphcoord),/data,/to_dev))([0,1],*)

index=where((NewPlanecoord(0,*) lt 1000) and (NewPlanecoord(1,*) lt 1000))

Flag=median_value-z

		if index(0) ge 0 then begin

z(NewPlanecoord(0,index),NewPlanecoord(1,index))=	$
	Data(Planecoord(0,index),Planecoord(1,index))

Flag(NewPlanecoord(0,index),NewPlanecoord(1,index))=255b

if keyword_set(outside) then Disk(NewPlanecoord(0,index),NewPlanecoord(1,index))=0
		endif

if keyword_set(show) then tvscl, z

index=0
NewPlanecoord=0
Planecoord=0

index=where(flag eq 0)

Flag=0
	if index(0) ge 0 then begin

Planecoord=(subs1to2(temporary(index),dimension=[Sz(1),Sz(2)]))

Sphcoord=(convert_coord(Planecoord,/dev,/to_data))([0,1],*)

Sphcoord(0,*)=Sphcoord(0,*) mod 360 < Latmax > Latmin
Sphcoord(1,*)=Sphcoord(1,*) mod 360 < Lonmax > Lonmin

if n_elements(LATITUDE) le 0 then $
	Sphcoord(0,*)=Sphcoord(0,*)+difrot(dt,Sphcoord(1,*),/days,/degree) $
		> Latmin < Latmax else $
	Sphcoord(0,*)=Sphcoord(0,*)+difrot(dt,Sphcoord(1,*),/days,/degree, lat=LATITUDE) $
		> Latmin < Latmax

NewPlanecoord=(convert_coord(temporary(Sphcoord),/data,/to_dev))([0,1],*)

z(Planecoord(0,*),Planecoord(1,*))=Data(NewPlanecoord(0,*),NewPlanecoord(1,*))

if keyword_set(show) then tvscl, z

NewPlanecoord=0
Planecoord=0

	endif

meridian=indgen(Latmax-Latmin+1)+Latmin

	if dt gt 0 then begin

if n_elements(LATITUDE) le 0 then $
;Edge=(convert_coord(-90 > Lonmin < Lonmax-difrot(dt,0,/days,/degree) mod 360,	$
;	meridian,/data,/to_dev))([0,1],*) else $
;Edge=(convert_coord(-90 > Lonmin < Lonmax-difrot(dt,0,/days,/degree, lat=LATITUDE) mod 360,	$
;	meridian,/data,/to_dev))([0,1],*)

range=[Nm/2, n_elements(xcir)-1]

;empty_field=polyfillv(	$
;	[xcir(range(0):range(1)),transpose(Edge(0,*))],		$
;	[ycir(range(0):range(1)),transpose(Edge(1,*))],		Sz(1),Sz(2))

	endif else begin

if n_elements(LATITUDE) le 0 then $
;Edge=(convert_coord(90 > Lonmin < Lonmax-difrot(dt,0,/days,/degree),	$
;	meridian,/data,/to_dev))([0,1],*) else $
;Edge=(convert_coord(90 > Lonmin < Lonmax-difrot(dt,0,/days,/degree, lat=LATITUDE),	$
;	meridian,/data,/to_dev))([0,1],*)

range=[0, Nm/2-1]

;empty_field=polyfillv(	$
;	[transpose(Edge(0,*)),reverse(xcir(range(0):range(1)),1)],	$
;	[transpose(Edge(1,*)),reverse(ycir(range(0):range(1)),1)],	Sz(1),Sz(2))

	endelse

;if empty_field(0) ge 0 then z(empty_field)=missing

;z(Edge(0,*), Edge(1,*))=missing
if keyword_set(show) then tvscl, z

if keyword_set(show) then empty

if not keyword_set(show) then wdelete,N_win

	if strmid(!version.release,0,1) ge 5 then begin

z1=z*0
z1(Disk_index)=z(Disk_index)
Disk_index=0
z=temporary(z1)

if keyword_set(show) then tvscl, z

	endif

if keyword_set(outside) then z=z+Data*Disk

return, z

		end